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ABSTRACT 


Spatial  competition  is  a  phenomenon  which  exists  in  numerous 
situations  throughout  nature.  Examples  of  objects  in  competition 
for  space  can  be  cited  for  biological,  economic,  social,  and 
environmental  processes.  In  many  of  these  processes,  the 
competition  is  for  resources  distributed  in  space,  rather  than 
simply  the  space  itself.  When  two  or  more  objects  compete  for  the 
same  space,  it  is  possible  to  model  the  process  such  that  we  may 
objectively  define  the  resulting  space  occupied  by  the  objects. 

In  this  paper,  I  develop  a  model  for  the  competitive  process  in 
order  to  define  the  region  of  space  occupied.  The  model  is 
incorporated  into  a  methodology  which  utilizes  a  coordinate 
transformation  function  to  determine  the  space  the  objects  occupy. 
An  application  of  the  methodology  to  groundwater  protection  is 
included  in  which  water  wells  are  the  objects  and  groundwater  is 
the  resource  being  competed  for.  The  study  was  conducted  using  a 
geographic  information  system  to  transform  spatial  coordinates  and 

generate  figures  defining  the  occupied  space.  101  pp. _ 
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CHAPTER  ONE 


INTRODUCTION 

Water  is  viewed  by  most  people  as  one  of  nature's  most 
valuable  resources.  Without  it  life  cannot  be  sustained. 
Supplies  of  water  are  often  mistakenly  perceived  as 
virtually  unlimited,  especially  in  relation  to  one's 
personal  needs.  It  is  crucial  to  note  that  the  availability 
of  an  adequate  water  supply  is  dependent  on  the  resource ' s 
distribution  through  space  and  time.  Even  when  available,  it 
can  be  made  unusable  if  polluted. 

Surface  water  and  ground  water  are  the  major  sources 
of  public  water  supply  in  the  United  States.  Both  sources 
are  susceptible  to  pollution.  Surface  supplies  are  readily 
contaminated  due  to  their  unconfined  nature  (lack  of 
barriers  between  the  supply  and  the  environment).  Ground 
water,  because  of  its  confined  character,  is  polluted 
through  more  subtle  processes.  Contaminated  supply  wells 
must  be  condemned  if  the  resource  is  determined  unusable. 

The  polluted  groundwater  can  even  contaminate  surface  water 
bodies  fed  by  the  associated  aquifer  (groundwater 
reservoir).  Due  to  the  relatively  slow  movement  of 
ground  water,  aquifer  contamination  may  not  be  detected 
until  years  after  a  pollutant  release. 
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Groundwater  Protection 

^Taylor  (1977)  suggests  that  source  protection  is  the 
"first  line  of  defense"  in  assuring  drinking  water  quality. 
The  intent  of  protective  measures  is  to  reduce  or  eliminate 
human  activities  which  utilize  substances  that  pose  a 
potential  contamination  threat  to  the  water  source.  — - 
Criteria  for  defining  such  substances  are  as  follows  (van 
Hageningh  and  van  Dui jvenbooden,  1979): 

dangerous  substances  are  substances  that 
are  harmful  to  health 

(-  impair  taste,  aroma,  or  color  of  water 

-  lead  to  lack  of  subsoil  aeration 

-  increase  the  corrosiveness  or  hardness  of  the 

water  _ 

-  affect  the  water  temperature 

Water  pollution  can/be  classified  by  contaminant  type. 
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biological  or  chemical.  >The  primary  aim  of  the  water  supply 
industry  has  been  to  ensure  that  the  consumer  is  provided 
with  water  which  is  free  from  pathogenic  microorganisms . - 
Yet,  the  single  most  prominent  cause  of  waterborne  disease 
in  the  U.S.  is  the  presence  of  pathogens  due  to  inadequate 
disinfecting  procedures  (Toft,  1985).  ^In  spite  of  these 
incidences  of  bacteriological  contamination,  the  major 
public  concern  with  drinking  water  has  been  possible 
contribution  to  cancer  risks  due  to  chemical  pollutants.  JJl  ~ 
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To  date,  over  8,000  wells  have  been  closed  or 


substantially  affected  by  toxic  contaminants  (Bachman,  1986) 


with  more  than  200  chemical  substances  detected  in  the 


groundwater  (Office  of  Technology  Assessment,  1985).  This 


number  is  expected  to  increase  dramatically.  The  expected 


increase  is  due  to  potential  exposure  of  ground  water  to  the 


more  than  60,000  known  chemicals  involved  in  manufacturing 


(Shackelford  and  Cline,  1986).  Of  this  huge  number  of 


potential  chemical  contaminants,  only  250  are  regulated 


under  the  Resource  Conservation  and  Recovery  Act,  and  even 


these  few  have  no  proven  methods  or  standards  for  quality 


assurance  or  quality  control  (Dowd,  1985). 


A  more  reasonable  and  flexible  approach  to  protection 


measures  is  to  list  those  activities  which  risk  release  of 


dangerous  substances.  Groundwater  contamination  processes 


can  be  tuvided  into  two  categories  based  on  the  activity 


which  resulted  in  the  pollution.  Those  activities  which  can 


be  located  precisely  are  known  as  point  source 


contamination.  Surface  impoundments,  landfills,  spills,  and 


leaks  are  examples  of  point  polluting  activities.  Diffuse 


activities,  such  as  agricultural  pesticide  and  fertilizer 


application,  and  residential  septic  tank  seepage  spread 


contaminants  over  large  areas.  These  activities  are  known 


as  nonpoint  pollution  sources.  Both  catagories  of 


contaminating  activities  are  equally  hazardous  to 
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groundwater  supplies,  but  the  latter  is  much  more  difficult 
to  isolate  and  identify. 

Florida's  Groundwater  Problems 

The  State  of  Florida,  which  has  more  ground  water  than 
any  other  state  (McGuinness,  1963),  supplies  88  percent  of 
the  public  water  consuming  households  from  groundwater 
sources  (Heath  and  Conover,  1981).  Although  the  aquifers 
containing  Florida's  groundwater  supply  extend  to  relatively 
great  depths,  the  shallow  zones  tapped  for  public  supply  are 
vulnerable  to  contamination.  The  overlying  materials  don't 
impede  the  movement  of  surface  contaminants  into  the 
groundwater  supply.  Considering  these  factors  and  the 
state's  intense  growth  and  development  trends,  one  can 
easily  see  the  critical  importance  of  protecting  the  quality 
of  Florida's  groundwater  resources.  The  increasing  demands 
on  a  vulnerable  groundwater  supply  necessitate  regulatory 
measures  to  assure  the  public  an  adequate  quantity  of 
potable  water. 

Recognizing  the  importance  of  regulatory  protection 
for  the  state's  groundwater  resources,  the  Florida 
Department  of  Environmental  Regulation  ( FDER )  developed 
specifications  for  delineating  "protection  zones"  around 
vulnerable  public  supply  wells  (Vecchioli  et  al . ,  1987). 
Within  these  protection  zones,  certain  activities  which 
entail  contamination  risks  to  the  ground  water  are 
prohibited  or  restricted  by  Florida  statutes.  As  a 
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consequence  of  the  legislation,  the  state  is  required  to  map 
the  protection  zones  around  existing  public  supply  wells. 

The  FDER  solicited  the  assistance  of  the  U.S.  Geological 
Survey  Water  Resources  Division  to  define  and  map  the  zones. 

Protection  Zone  Definition 

Fried  and  Zampettl  (1979,  p.292)  defined  a  groundwater 
protection  zone  as  "a  geographically  circumscribed  area 
where  certain  activities  are  banned,  controlled,  or  subject 
to  authorization,  and  the  use  and  ownership  of  which  are 
subject  to  certain  obligations  in  order  to  safeguard  the 
underground  waters  against  the  risk  of  pollution."  These 
zones  are  restrictive  areas  situated  about  supply  wells 
(Figure  1)  with  the  intent  of  limiting  surface  land  use 
activities  to  preserve  the  quality  of  the  underlying  water 
reserves  contained  in  the  aquifer.  The  protection  of  an 
aquifer  is  a  popular  notion,  but  it  is  a  much  more  complex 
problem  than  is  initially  apparent.  Difficulties  arise  when 
an  attempt  is  made  to  precisely  define  what  surface  area 
affects  the  groundwater  supply.  The  major  problem  is  a  lack 
of  scientific  understanding  of  the  movement  of  water  and  the 
attenuation  of  contaminants  in  the  ground.  Site  specific 
studies  to  determine  the  geohydrologic  characteristics  of  an 
area  are  few  in  number.  Research  of  this  nature  is 
expensive  and  frequently  beyond  the  means  of  local  agencies. 
This  leads  to  a  general  problem  often  confronted  when  one 
wishes  to  Impose  environmental  safeguards:  the 
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GROUNDWATER  PROTECTION  ZONE 


WELL 


CALCULATED  5  YEAR  AQUIFER  VOLUME 


Figure  1:  Cross-Section  Of  A  Groundwater  Protection  Zone 
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cause-and-ef f ect  relationship  may  not  be  easily  defined.  In 
the  case  of  ground  water,  this  i«  especially  true  since  the 
water,  soil,  and  pollutant  are  all  hidden  from  view. 
Therefore,  as  Goldrosen  (1977,  p.57)  has  stated,  standards 
and  rules  for  groundwater  protection  are  "still  a  matter  of 
guesswork . " 

When  attempting  to  regulate  activities  near  a 
withdrawal  point,  the  uncertainties  associated  with  the 
groundwater  pollution  process  tend  to  complicate  matters. 
Initially,  one  might  assume  that  the  protective  measures 
should  cover  the  entire  catchment  area.  This  seems 
reasonable  since  a  polluting  activity  located  within  an 
aquifer  intake  area  could  potentially  contaminate  any  well 
drawing  from  that  aquifer.  However,  these  aquifers  can 
cover  vast  areas.  Protective  measures  require  that 
limitations  be  placed  on  land  use,  resulting  in  reductions 
in  land  value.  The  individual  landowners  whose  property 
values  may  be  diminished  will  certainly  be  reluctant  to 
support  such  legislation.  Despite  the  fact  that  the 
principle  of  protecting  the  entire  intake  area  is 
theoretically  correct,  deviations  from  it  are  necessary  for 
practicality.  The  protection  zones  need  to  be  as  large  as 
required  for  health  and  safety  purposes,  but  as  small  as 
possible  for  socio-economic  reasons.  Only  when  a  fair 
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compromise  is  reached  between  these  two  contradictory 
factors  can  protective  legislation  be  enacted. 


CHAPTER  TWO 


PROBLEM  STATEMENT 

Groundwater  protection  zone  delineation  is  a  problem 
which  can  easily  be  recognized  at  the  practical  level.  To 
get  a  greater  understanding  of  the  problem,  one  must 
consider  it  at  the  conceptual  level  as  well.  Additional 
insight  and  scientific  utility  are  derived  from  a 
fundamental  theoretical  approach.  For  this  reason,  this 
investigation  addresses  both  the  conceptual  and  the  applied 
aspects  of  this  spatial  delineation  problem. 

Description  of  the  problem  in  an  abstract  manner  is 
somewhat  awkward  and  undoubtedly  confusing  without  mention 
of  its  practical  application.  The  perception  of  the  spatial 
concepts  involved  becomes  more  concrete  when  one  first 
understands  the  circumstances  which  initially  pose  the 
problem.  The  protection  zone  delineation  problem  is 
described  first.  This  should  facilitate  greater 
comprehension  of  the  description  and  implications  of  the 
fundamental  theoretical  concepts. 

Practical  Problem 

FDER  rule  docket  number  85-23R  specifies  the 
requirements  for  delineation  cf  protection  zones  in  the 
State  of  Florida  (Vecchioli  et  al.,  1987).  The  size  of  the 
zone  for  an  individual  well  is  calculated  using  a  radial 
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plug-flow  displacement  model  based  on  a  five-year 
permltted-rate  withdrawal.  The  FDER  formula  used  to 
determine  the  zone  size  is: 

r  =  [  Q*t/7T*h*n  ]o.5 

where 

r  ■  radius  of  the  individual  protection  zone 
t  -  time,  5  years  or  1825  days 

n  =  effective  porosity  of  the  associated  aquifer 


a 


i 


TT  =*  3.14 

Q  =  permitted  withdrawal  rate  (ft3/day)  ,and 
h  »  thickness  of  the  aquifer  interval  penetrated  (ft) 
The  result  of  the  calculation  is  the  radius  of  a  circle. 

When  circumscribed  about  the  supply  well,  this  circle 
defines  the  zone  of  protection  for  that  well. 

In  some  instances,  the  protection  zones  of  neighboring 
wells  overlap.  This  is  very  evident  in  areas  where  a  number 
of  wells  are  concentrated  in  space,  known  as  wellfields. 

The  formula  for  determining  the  size  of  the  protection  zones 
is  based  on  a  volumetric  displacement  concept.  The  volume 
of  the  water  to  be  protected  is  proportional  to  the  area  of 
the  protection  zone.  Thus,  a  condition  of  overlapping  zones 
results  in  underestimation  of  the  area  requiring  protective 
measures.  The  amount  of  area  by  which  the  overlapping  zones 
are  deficient  for  groundwater  protection  is  equivalent  to 
the  area  of  overlap. 
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The  FDER  regulation  stipulates  that  for  wells  whose 
zones  overlap: 

"...the  area  on  the  surface  overlying  the  aquifer  equal  to 
the  sum  of  the  areas  of  the  five  year  zones  of  protection  of 
the  individual  wells,  shall  be  used  to  define  the  area  which 
encircles  the  perimeter  of  the  wellfield.  In  cases  where  a 
zone  of  protection  of  a  single  well  protrudes  beyond  the 
calculated  perimeter  or  when  the  configuration  of  the 
wellfield  is  irregular,  the  perimeter  will  be  shaped  to 
accommodate  the  configuration.  The  surface  area  encircling 
the  perimeter  of  the  wellfield  shall  not  exceed  the  total 
surface  areas  of  the  overlapping  zones  of  protection  for 
individual  wells"  (Vecchioli  et  al . ,  1987,  p.60). 

The  U.S.  Geological  Survey  procedure  for  delineating 
the  zones  for  the  overlap  condition  requires  extending  the 
area  defined  by  the  union  of  the  overlapping  individual 
protection  zones  (Figure  2)  to  an  enlarged  area  which 
approximates  the  shape  of  the  well  configuration.  Through 
successive  approximations  the  area  is  expanded  until  it  is 
within  five  percent  of  the  summed  areas  of  the  individual 
zones.  If  upon  expansion  of  the  overlapping  zones,  the 
boundary  of  the  enlarged  area  contacts  a  previously 
non-overlapping  protection  zone,  the  contacted  zone  is  added 
to  the  union  and  the  total  areal  requirement  recomputed. 
Successive  expansions  produce  a  figure  of  correct  area  whose 
shape  is  representative  of  the  well  configuration  (Figure 
3).  This  expanded  region  is  delineated  on  its  respective 
quadrangle  sheet.  The  resultant  area  mapped  is  designated  a 
composite  protection  zone. 

As  is  evident  from  the  above  description,  the  U.S. 
Geological  Survey's  manual  cartographic  procedure  for 


UNION  OF  OVERLAPPING 
PROTECTION  ZONES 


Figure  2:  Union  Of  Overlapping  Protection  Zones 
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EXPANDED  CONFIGURATION  FOR 
COMPOSITE  PROTECTION  ZONE 


Figure  3:  Expanded  Figure  For  Composite  Protection  Zone 
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defining  composite  protection  zones  is  subjective.  The 
process  relies  on  an  individual's  judgement  to  delineate  a 
region  shaped  like  the  well  field  configuration.  Different 
analysts  could  arrive  at  dissimilar  composite  protection 
zones.  Since  the  purpose  of  protection  zone  mapping  is  to 
define  land  area  upon  which  land  use  restrictions  will  be 
applied,  this  subjectivity  poses  problems  to  the  legal 
aspects  of  enforcement.  Owners  of  use-restricted  land 
could  certainly  challenge  the  "arbitrary"  method  by  which 
the  areas  were  defined. 


Conceptual  Aspects 

At  the  conceptual  level ,  the  problem  equates  to  the 
modeling  of  a  spatial  delineation  process.  This,  like  most 
other  spatial  processes,  is  characterized  using  a  general 
model  framework  employing  points,  lines,  and/or  areas. 
Patterns  of  these  three  feature  types  can  be  generated  by 
the  model  according  to  various  rules.  These  rules  represent 
assumptions  inherent  to  the  spatial  process  and  thus  define 
the  spatial  process  model. 

The  patterns  of  concern  here  are  area  patterns 
resulting  from  the  spread  of  influence  away  from  a  set  of 
Initial  point  sources.  This  type  of  spatial  process  model 
is  commonly  termed  a  growing  space  model.  The  growing  space 
model  employs  two  basic  concepts  of  importance,  competition 


and  space. 


Cormack  (1979,  p.152)  defines  competition  as  "any 
process-physical  inhibition,  direct  interference,  or 
overlapping  of  niches-which  has  a  negative  influence  on  an 
individual's  location,  viability,  or  strength."  He  further 
describes  space  as  "the  two-dimensional  physical  area 
required  by  an  organism  for  development."  The  natural 
representation  of  space  occupied  by  an  organism  is  assumed 
to  be  a  circle  centered  on  the  organism.  If  the  environment 
is  uniform,  the  space  occupied  can  also  represent  the 
resources  available. 

Obviously,  the  terms  used  by  Cormack  are  indicative  of 
biological  processes,  but  the  concepts  involved  are 
evidently  appropriate  to  the  present  modeling  problem. 
Competition  for  space/resources  can  result  in  two  distinct 
Interactions.  First,  the  "growing"  circles  continue  to  grow 
after  contact  with  neighboring  circles.  Thus,  overlap 
occurs.  The  amount  of  overlap  reflects  the  competition 
pressure  imposed  by  one  organism  on  another.  Alternatively, 
the  growing  circles  are  not  allowed  to  overlap,  but 
distorted  into  polygons  which  represent  the  space  available 
to  the  organism  for  development. 

Getis  and  Boots  (1978)  describe  the  growth  model  in 
more  general  terms  and  refer  to  spatial  processes  of  this 
type  as  area  generating  point  processes.  Those  processes 
which  prohibit  overlap  generate  patterns  denoted  as 
contiguous.  This  type  of  process  leads  to  a  pattern  which 


is  space-exhaustive  (Figure  4).  The  process  type  which 
permits  overlap  generates  patterns  designated  as  overlap 
patterns,  and  is  not  space-exhaustive. 

The  present  study  is  concerned  with  a  process  which  is 
a  hybrid  of  the  two  types  described  previously.  Initially, 
overlap  is  allowed  as  a  means  of  identifying  the  amount  of 
overlap  present.  After  the  size  of  the  overlap  area  is 
determined,  the  alternate  condition  which  prohibits  overlap 
is  imposed.  An  area  equal  to  the  overlap  area  then  must  be 
distributed  about  the  overlap  pattern  generated.  Modeling 
this  allocation  process  is  the  primary  focus  of  this  study. 
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Figure  4:  Growth  Patterns 
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CHAPTER  THREE 


STUDY  OBJECTIVES 

The  intent  of  this  study  is  to  develop  a  model  for  a 
process  of  spatial  delineation  based  on  the  concepts  of 
competition  and  space  described  earlier.  From  this  model, a 
methodology  is  created  for  composite  protection  zone 
delineation  which  attempts  to  minimize  scientific  and 
cartographic  subjectivity.  As  a  means  of  satisfying  this 
intent,  the  basic  criteria  defining  the  goals  of  the  model 
must  be  stated.  It  seems  most  useful  to  outline  the  goals 
of  the  modeling  process  in  the  context  of  its  application. 

First ,  the  model  must  meet  the  two  requirements 
explicitly  stated  in  the  FDER  rule.  The  composite 
protection  zone  must  totally  enclose  all  the  circular  zones 
of  the  individual  wells.  The  area  of  the  composite  zone 
should  be  equal  to  the  sum  of  the  areas  of  the  individual 
zones.  Both  requirements  are  met  by  the  U.S.  Geological 
Survey  procedure,  although  the  total  area  is  to  a  lesser 
degree  (5  percent  tolerance  level). 

At  the  heart  of  the  matter  lies  the  objectivity  of  the 
process.  The  fundamental  reason  for  the  study  is  to  develop 
a  formal  approach  to  the  delineation  problem.  Such  a 
procedure  should  yield  a  unique  composite  zone  which  can  be 
reproduced  during  independent  replication  of  the  process. 
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Only  when  independent  attempts  result  in  identical 
protection  zones,  will  the  objectivity  of  the  methodology  be 
substantiated . 

The  underlying  problem  which  introduces  subjectivity 
into  the  manual  procedure  is  that  decisions  concerning  the 
extension  of  the  protection  zone  boundaries  are  based  on 
visual  inspection  of  the  overlapping  well  configuration. 

The  key  to  solving  this  problem  lies  in  the  ability  to 
define  strict  rules  which  dictate  precisely  how  the  overlap 
area  of  two  individual  protection  zones  is  distributed  when 
extending  zone  boundaries. 

An  attempt  was  made  to  quantify  the  decision  process 
utilized  during  manual  delineation.  U.S.  Geological  Survey 
personnel  were  queried  concerning  cognitive  guidelines  which 
governed  the  distribution  of  the  overlap  while  generating 
composite  zones.  Though  difficult  to  define  and  more 
difficult  to  quantify,  the  discussion  with  USGS  personnel 
resulted  in  one  general  conclusion:  the  overlapping  zones 
were  extended  based  on  proximity  to  the  center  of  the 
overlap  area.  The  zone  boundary  was  expanded  more  near  the 
center  of  the  overlap  and  to  a  lesser  degree  at  more  distant 
points.  This  is  an  inverse  relationship  between  distance 
from  overlap  center  and  boundary  extension.  An  inverse 
function  based  on  distance  from  overlap  center  seemed  to  be 
the  most  plausible  and  reasonable  approach  to  the  problem, 
according  to  Geological  Survey  hydrologists.  Central  to 


this  study,  the  premise  of  an  inverse  distance  relation  is 
assumed  to  be  accurate,  based  on  U.S.  Geological  Survey 
expertise  of  groundwater  movement. 

The  study  results  in  a  methodology  developed  from  a 
model  incorporating  this  inverse  distance  relationship. 
Testing  requires  the  generation  of  composite  protection 
zones  in  a  format  similar  to  those  presented  in  the 
U . S . Geological  Survey  Report  88-4051  (Vecchioli  et  al . , 

1987) .  An  analysis  and  evaluation  is  also  conducted  to 
establish  the  acceptability  of  the  approach  for  application 
to  groundwater  protection  as  well  as  other  spatial 
competition  problems. 

In  summary,  the  objectives  of  this  study  are  to 
develop  and  evaluate  a  methodology  for  generating  composite 
protection  zones  based  on  a  spatial  distribution  model.  The 
proposed  method  is  subject  to  the  following  constraints: 

-  the  composite  protection  zone  must  totally  enclose 

all  individual  protection  zones  of  which  it  is 
comprised 

-  the  area  of  the  composite  zone  must  equal  the  sum  of 

the  areas  of  the  individual  protection  zones  of 
which  it  is  comprised 

-  the  method  must  be  formal  and  yield  a  unique 

composite  protection  zone  of  reproducible  nature 


-  the  central  premise  concerning  overlap  distribution 
(boundary  extension)  will  incorporate  an  inverse 
relationship  between  distance  from  overlap  center  to 
boundary  and  boundary  extension. 

Utilizing  a  delineation  methodology  employing  these 
constraints,  objective  composite  protection  zone  delineation 
should  be  possible.  In  addition,  the  model  may  be 
applicable  to  similar  space  allocation  problems. 


CHAPTER  POUR 


RELATED  LITERATURE 


A  review  of  pertinent  literature  is  included  to  address 
the  present  state  of  knowledge  in  collateral  subject  areas. 
The  review  includes  references  to  groundwater  protection  as 
well  as  to  relevant  spatial  modeling.  Previous  work  which 
lends  to  the  study  effort  are  utilized  to  the  greatest 
possible  advantage  toward  problem  solution. 

Groundwater  Issues 

The  increasing  level  of  concern  for  the  protection  of 
ground  water  in  the  United  States  is  apparent  from  the 
amount  of  legislation  enacted  to  protect  the  resource.  In 
addition  to  Florida,  the  states  of  Arizona,  California, 
Colorado,  Connecticut,  Kansas,  Massachusetts,  New  Jersey, 

New  York,  and  Wisconsin  have  passed  comprehensive 
groundwater  protection  statutes  (Mosher,  1987).  Still,  the 
Florida  program  to  protect  zones  around  wells  with  land  use 
controls  is  a  unique  approach. 

The  control  of  land  use  for  groundwater  protection  is 
not  without  controversy.  The  lack  of  consistent, 
comprehensive  data  has  made  it  difficult  to  establish  direct 
relationships  between  human  activities  and  groundwater 
contamination  (Lee  and  Nielsen,  1987).  The  characteristics 
of  groundwater  movement  are  poorly  understood  by  the  public. 


Frequently,  the  individual  owning  land  with  restrictions 
applied  to  its  use  will  react  negatively  to  such  regulation. 
Libby  and  Kovan  (1987)  indicate  that  the  lack  of  general 
awareness  about  connections  between  human  activities  and 
groundwater  quality  are  responsible  for  such  attitudes. 

Land  use  controls  are  an  essential  means  of 
groundwater  protection,  since  discharges  from  nonpoint 
sources  are  not  amenable  to  technological  controls  (Pye  et 
al.,  1983).  Therefore,  these  sources  of  contamination  must 
be  controlled  locally.  The  primary  need  in  local  land  use 
planning  is  the  definition  of  areas  of  influence  for 
operating  wells  (Libby  and  Kovan,  1987).  Once  the  areas  are 
identified,  land  use  patterns  in  a  community  can  be  guided 
in  ways  which  enhance  groundwater  protection. 

Land  use  controls  are  a  means  of  segregating 
contaminative  activities  from  the  water  supply.  Pollution 
sources  should  be  as  remote  as  possible  from  the  groundwater 
reservoir  and,  hence,  withdrawal  and  catchment  areas.  This 
gives  rise  to  protective  areas  or  "zones"  enclosing 
abstraction  points  and  catchment  areas. 

Several  zones  are  usually  established  (Fried  and 
Zampetti,  1979).  An  inner  zone  of  a  relatively  small  area 
around  the  well  is  defined.  In  this  inner  zone  all 
activities  with  contamination  potential  are  banned. 

Normally  the  water  authority  is  required  to  own  land  within 
this  zone.  An  intermediate  zone  is  specified  which  is 


ry  ■>  'J*  r'jrrjrTjr 


23 

designed  to  protect  against  bacteriological  contamination. 
Within  this  area,  activities  which  might  introduce 
biological  toxins  to  the  groundwater,  such  as  cattle  grazing 
or  septic  leaching,  are  prohibited.  An  outer  zone  is 
designated  to  cover  the  entire  catchment  area.  Within  this 
outer  zone  potential  chemical  polluters  are  banned. 
Controlling  this  larger  area  allows  additional  time  for 
intervention,  dilution,  and  breakdown  of  chemical 
contaminants . 

The  primary  method  for  determining  the  size  of 
protection  zones  is  one  based  on  delay  times.  Delay  time 
refers  to  the  time  it  takes  a  water  droplet  located  at  a 
given  distance  from  a  well  to  travel  to  the  same  well  under 
forces  of  water  withdrawal.  Calculating  zone  size  utilizing 
delay  time,  one  can  take  into  account  the  geohydrology  of 
the  subsoil,  and  size,  form,  and  capacity  of  the  well  (van 
Wageningh  and  van  Dui jvenbooden,  1979). 

Langeweg  (1979)  discusses  a  three-zone  system  used  as  a 
guideline  in  the  Federal  Republic  of  Germany.  It  consists 
of  zones  numbered  I, II,  and  III,  varying  in  size  from  10 
meters,  50-60  days  delay  time,  and  the  entire  intake  area, 
respectively.  The  Netherlands  also  uses  a  three- 
protection-zone  system  for  groundwater  wells  (van 
Di jk-Looi jaard  and  de  Kruijf,  1985).  An  inner  60-day  zone 
is  required  to  prevent  bacteriological/chemical 
contamination  and  bans  all  activity  within  its  boundary.  A 
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10-year  intermediate  zone  is  specified  whose  delay  time 
chosen  to  be  adequate  to  control  pollution  after  an  accident 
or  if  necessary,  replace  the  well.  The  outermost  zone  of 
25-year  delay  time  is  based  on  the  time  required  to  replace 
a  well  and  treatment  plant  without  financial  backlash 
(adequate  write-off  period).  The  restrictiveness  of  the 
Dutch  regulation  decreases  from  inner  to  outer  zones,  with 
the  10  and  25-year  zones  addressing  potential  chemical 
polluting  activities  only. 

In  the  United  States,  legislation  concerning 
groundwater  protection  zones  has  been  much  less  ambitious. 
Many  state  regulations  require  only  that  water  supply  wells 
be  located  100  feet  or  more  from  potential  sources  of 
pollution  (USEPA,  1985).  The  leniency  of  the  environmental 
regulations  can  in  part  be  traced  to  the  two  major  drawbacks 
to  groundwater  classification  (zone  delineation)  posed  by 
Duda  and  Johnson  (1987).  First,  analyses  of  this  nature 
require  extensive  hydrogeologic  studies  which  may  upon 
initial  reflection  seem  too  costly.  Secondly,  there  exist  a 
number  of  technical  uncertainties  in  establishing  boundaries 
of  the  classified  units.  Both  shortcomings  must  be  resolved 
if  the  states  are  to  enact  regulatory  controls  which  will 
satisfy  their  objectives. 

Admittedly,  the  method  the  FDER  has  chosen  to  calculate 
and  delineate  protection  zones  in  Florida  is  not  the  most 
sophisticated  available.  Its  limitations  arise  from  the 
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lack  of  consideration  of  the  natural  slope  of  the 
potentiometric  surface,  recharge,  natural  boundaries, 
porosity  variations,  and  incomplete  knowledge  of  aquifer 
thickness.  Data  availability,  combined  with  resource 
constraints  for  calculating  and  mapping  the  zones,  dictated 
the  method  chosen.  Nevertheless,  the  technique  selected 
provides  a  major  step  toward  ensuring  Florida  of  an  adequate 
drinking  water  supply. 

Merchant  et  al.  (1987)  consider  their  approach  to 
capture  (protection)  zone  approximation  to  be  a  great 
improvement  over  the  circular  zone  method.  Taking  into 
account  regional  flow,  their  method  produces  elliptical 
zones.  By  generalizing  the  aquifer  to  a  homogeneous 
uniform-flow  system,  elliptical-shaped  zones  are  generated 
whose  eccentricity  increases  as  the  regional  flow  rate  at 
the  well  increases.  In  order  to  use  this  method,  more 
site-specific  information  would  be  required  than  for  the 
FDER  approach,  including  knowledge  of  the  local  flow 
direction  and  hydraulic  gradient.  Presently,  this  is  a 
prohibitive  factor  when  one  considers  zone  determination  for 
the  large  number  of  wells  involved  in  Florida. 

Spatial  Concepts 

The  problem  of  protection  zone  determination  has  been 
addressed  only  in  the  context  of  groundwater  application. 

At  a  more  fundamental  level,  it  can  be  viewed  as  a 
two-dimensional  spatial  problem.  Essentially,  two  "objects" 


located  at  specific  sites  are  in  competition  for  "resources" 
distributed  in  space.  Whether  the  objects  are  bacteria, 
plants,  or  wells,  and  the  resources  are  nutrients,  light,  or 
water,  the  goal  is  to  define  the  region  in  space  the 
competing  objects  require  to  fulfill  their  respective 
resource  needs.  For  the  purpose  of  this  study,  the  solution 
is  subject  to  the  constraints  previously  listed.  The 
resulting  solution  then  provides  the  shape  of  the  resource 
space  required  by  the  objects. 

Johnson  and  Mehl  (1939)  studied  surficial  processes  of 
the  formation  and  growth  of  metallic  nuclei.  They  described 
the  crystallization  of  the  nuclei  as  a  radial  growth 
process.  In  addition  to  their  work  concerning  changes  in 
growth  rate  of  the  crystal  structure  through  time,  they 
discuss  what  they  call  the  Impingement  factor.  This  factor 
is  so  named  because  it  describes  the  change  in  growth  rate 
produced  when  a  growing  crystal  contacts  or  impinges  on 
neighboring  crystals.  Sometimes  referred  to  as  contact 
inhibition,  this  factor  is  calculated  as  the  proportion  of 
the  total  surface  area  outside  the  crystal  considered  which 
is  not  covered  by  other  crystals.  It  is  used  to  determine 
the  fraction  of  a  plane  surface  covered  by  growth  with 
Impingement  through  time. 

In  a  similar  study  by  Evans  (1945),  the  covering  of  a 
surface  by  metallic  films  or  layers  of  corrosion  is 
examined.  He  describes  the  film  or  layer  spreading  out  as 


expanding  circles  from  predetermined  nuclei  distributed 
randomly  by  a  Poisson  process.  He  also  derives  an  equation 
for  calculating  the  area  of  a  surface  not  covered  by  the 
circles  which  is  based  on  nuclei  density,  radial  growth 
velocity,  and  time.  The  equation  is  then  modified  to 
consider  nucleation  rate,  nucleation  number,  and  three- 
dimensional  growth.  Evans  further  expands  his  idea  to 
growths  which  avoid  each  other,  such  as  is  the  case  when  a 
substance  necessary  for  growth  becomes  exhausted  as  a  result 
of  growth.  Unfortunately,  he  concludes  discussion  of  this 
matter  suggesting  that  this  situation  "often  leads  to 
dendritic  structures"  (p.368). 

Studies  with  potential  implications  can  also  be  found 
among  the  more  conventional  spatial  analysis  approaches. 
Nearest  neighbor  methods  cover  a  variety  of  techniques  which 
use  distance  measures  to  identify  and  analyze  point 
patterns.  The  ultimate  goal  of  these  methods  is  to 
understand  the  process  that  generates  the  patterns.  By 
definition,  a  point  nearest  to  a  randomly  chosen  point  is 
its  nearest  neighbor.  The  distance  between  the  two  points 
is  the  nearest  neighbor  distance.  It  is  this  distance  that 
is  used  in  nearest  neighbor  studies. 

A  commonly  used  technique  is  the  Clark  and  Evans 
nearest  neighbor  procedure  (Ripley,  1981).  This  approach 
compares  the  average  nearest  neighbor  distance  with  the 
expected  average  distance  for  randomly  placed  points.  In 


order  to  determine  the  expected  distance,  one  must  know  the 
point  density  (number  of  points  per  unit  area) .  Normally, 
ascertaining  the  number  of  points  in  consideration  is  not  a 
significant  problem.  But  how  to  define  the  area  of  the 
study  region  can  be  quite  a  different  matter  when  dealing 
simply  with  a  set  of  points.  Numerous  suggestions  have  been 
presented  as  possible  approaches,  but  none  are  universally 
accepted.  The  area  boundary/edge  effect  problem,  as  it  is 
referred  to,  is  the  major  difficulty  encountered  with  all 
nearest  neighbor  techniques  (Upton  and  Fingleton,  1985). 

The  boundary  effect  associated  with  these  techniques  is  a 
problem  similar  to  the  areal  delineation  problem  of  this 
study.  However,  the  boundary  solutions  suggested  for 
nearest  neighbor  methods,  such  as  buffer  zones  or  toroidal 
surfaces  (Ripley,  1981),  are  inappropriate  for  the  present 
investigation. 

Simberloff  (1979)  modified  the  Clark  and  Evans 
technique  to  consider  the  data  as  a  set  of  circles  centered 
at  the  points  identified,  rather  than  simply  as  points.  He 
was  concerned  with  plant  inhibition  models  which  viewed 
plants  as  non-overlapping  discs.  His  approach  allows  for 
calculation  of  the  average  radius  of  the  discs  based  on 
nearest  neighbor  distances.  Upton  and  Fingleton  (1985) 
provide  an  example  of  Simberloff 's  modification  in  an 
application  to  the  spatial  distribution  of  beadlet  anemones 
on  the  face  of  a  boulder.  Using  the  modification,  they  were 


able  to  estimate  the  radius  of  a  disc  occupied  by  an 
anemone.  A  comparison  of  the  organism's  size  with  that  of 
the  estimated  disc  radius  revealed  interesting  colonization 
characteristics.  The  aggressive  little  beasts  occupied  a 
territory  over  12  times  their  own  size.  However 
interesting,  this  approach  is  inappropriate  for  this  study, 
since  the  results  are  based  solely  on  average  interpoint 
distance. 

The  majority  of  inhibition  models, such  as  Simberlof f 1 s , 
were  principally  concerned  with  ecological  applications. 
However,  regular  patterns  are  also  known  to  exist  in  human 
geography.  Most  of  these  seek  deterministic  rather  than 
stochastic  explanations  for  the  resulting  patterns. 

A  dominant  theme  in  studies  of  this  type  is 
Christaller's  centra 1  place  model.  His  model  concerns  the 
optimum  placement  of  supply  points  for  goods  and  services 
required  by  a  surrounding  population.  An  inhibitory 
mechanism  prescribes  regularity  in  the  locational  pattern  of 
central  places.  There  are  a  number  of  levels  of  central 
places  based  on  the  number  of  locations  required  to  supply 
the  demand  for  any  particular  good  or  service. 

Assuming  a  uniform  environment,  a  reasonable  consumer 
would  travel  the  minimum  distance  to  a  central  place  to 
acquire  goods  or  services  supplied  there.  To  achieve  the 
distance  minimization,  the  highest  level  centers  are  ideally 
arranged  in  a  triangular  lattice  pattern  with  lower  level 


centers  similarly  interwoven.  The  resulting  geometry 
suggests  that  each  central  place  is  located  at  the  center  of 
a  hexagonal  market  area.  Losch  (1954)  was  able  to  show 
mathematically  that  hexagons  provided  "the  best  possible 
packing"  of  trade  areas  in  a  region  for  the  mutual  advantage 
of  the  consumer  and  producer.  The  regular  subdivision  of 
the  region  into  a  hexagonal  tessellation  pattern  was  the 
ideal  structure  according  to  the  central  place  model. 

Another  tessellation  which  receives  widespread  use  in 
spatial  studies  is  known  as  a  Dirichlet  tessellation.  This 
tessellation  has  been  used  in  a  number  of  contexts  and  goes 
by  a  variety  of  names,  including  Voronoi  polygons,  Thiessen 
polygons,  the  cell  model,  and  the  S-mosaic  (Boots  and 
Murdoch.  1983).  The  polygons  which  comprise  the 
tessellation  are  formed  about  a  set  of  points  in  a  planar 
region.  Each  polygon  is  formed  about  a  single  point  and  has 
the  property  that  every  location  in  that  polygon  is  nearer 
to  that  point  than  any  other  point  in  the  region.  This 
results  in  subdivision  of  the  region  into  polygons  in  a 
space-exhaustive  manner  (Figure  5). 

Crain  (1972)  posed  that  stochastic  processes  in 
metallurgy,  cell  biology,  astrophysics,  and  geology  could 
all  be  modeled  using  the  Dirichlet  tessellation.  He 
described  an  example  concerning  two-dimensional  crystal 
growth  about  randomly  distributed  nuclei.  The  crystals  grow 
uniformly  over  the  plane  until  mutual  contact  prevents 
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further  growth.  The  process  results  in  coverage  of  the 
plane  with  non-overlapping  polygons  (Dirichlet  cells).  He 
cited  as  the  primary  interest  of  his  model  the  frequency 
distributions  of  such  polygon  parameters  as  the  number  of 
sides,  perimeter,  and  area.  Expected  values  were  also 
derived  for  each  of  these  parameters. 

The  basic  characteristics  of  the  polygons  formed  by  a 
Dirichlet  tessellation  are  given  by  Getis  and  Boots  (1978): 
(1)  all  polygons  have  straight  line  edges,  (2)  all  polygons 
are  convex  (i.e.,  it  is  possible  to  link  any  two  points 
within  the  polygon  without  traversing  another  cell),  and  (3) 
the  minimum  number  of  edges  is  three. 

If  Crain's  crystal  model  is  modified  to  incorporate 
different  rates  of  growth  between  individual  crystals  or  to 
allow  crystal  nuclei  to  be  "born"  at  different  times,  the 
resulting  pattern  is  described  as  a  curved  boundary 
tessellation  (Figure  6).  Such  a  process  of  crystal  growth 
was  described  previously  by  Johnson  and  Mehl  (1939).  The 
polygon  characteristics  for  this  model  differ  from  Dirichlet 
cells  in  that:  (1)  they  have  curvilinear  and  straight 
edges,  (2)  the  polygons  are  not  necessarily  convex,  (3)  the 
minimum  number  of  edges  is  two,  and  (4)  the  common  edge 
between  two  polygons  can  be  discontinuous. 

Other  than  the  above  characteristics,  very  little  is 
known  about  properties  of  models  of  this  type  (Upton  and 
Fingleton,  1985).  Recent  investigations  into  the  geographic 
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application  of  this  tessellation  have  yielded  somewhat 
inconclusive  results.  Getis  and  Boots  (1978,  p.147)  purport 
that  "because  of  its  potential  value  to  geographers,  it  is 
unfortunate  that  it  is  difficult  to  derive  analytically  the 
properties  of  patterns  produced  by  such  a  model". 

Geographical  research  involved  with  the  description  of 
shape  also  has  potential  implications  to  this  study.  An 
infinite  number  of  shapes  can  be  delineated  in  the  spatial 
competition  problem.  But  within  the  constraints  defined,  is 
there  a  unique  shape  which  delineates  the  protection  zone? 
Lee  and  Sallee  (1970)  suggest  that  a  finite  number  of 
parameters  is  not  sufficient  to  produce  a  unique  shape 
description.  Frolov  (1974)  indicates  that  these  parameters 
fail  to  adequately  account  for  characteristics  of  shape 
variation  which  he  denotes  as  compactness,  dissection,  and 
indentation.  The  parameters  referred  to  in  their  research 
are  single  value  types  such  as  area  and  perimeter.  Within 
the  constraints  employed  in  this  study,  there  is  a  single 
value  parameter  (area),  as  well  as  conditional  statements 
about  minimal  boundaries  and  overlap  distribution.  Since 
the  conditions  imposed  here  ar-  shape-definitive  in  nature, 
the  author  contends  that  unique  shapes  can  be  generated 
satisfying  all  constraints. 

In  reviewing  other  approaches  to  similar  spatial 
problems  for  potential  application  here,  the  works  by  Tobler 
and  Rushton  were  discovered.  Tobler  (1963)  described 


Hagerstrand' s  problem  of  depicting  migratory  exchange. 
Hagerstrand  suggested  that  migration  was  based  on  a  distance 
decay  function  with  decreasing  migratory  exchange  at 
increasing  distance  from  the  population  center.  To 
illustrate  this  he  used  a  logarithmic  cartogram  derived  from 
the  distance  function  resulting  in  radial  scale  distortion. 
Application  of  this  illustration  technique  to  the 
groundwater  problem  was  considered,  but  subsequently 
rejected  after  further  investigation.  It  was  discovered 
that  the  radial  scale  distortion  using  a  logarithmic 
function  extended  points  located  near  the  center  point  more 
than  those  further  from  the  center.  However,  those  points 
located  farthest  from  the  center  could  potentially  be  drawn 
in  toward  the  center,  a  fact  which  one  cannot  ignore 
considering  the  minimal  boundary  constraint.  However, 

Tobler  does  indicate  that  projections  similar  to 
Hagerstrand 1 s  can  be  obtained  by  using  suspected  functions 
of  radial  scale  distortion  for  any  distance  model  employed 
by  geographers. 

Rushton  (1972)  dealt  with  a  spatial  problem  as  a 
variation  of  the  central  place  model.  As  previously 
stated,  Christaller  assumed  environmental  conditions  to  be 
uniform.  Therefore,  he  suggested  that  ideally  located 
central  places  would  be  distributed  in  a  regular  pr  ttern. 


Rushton  argues  that  if  environmental  conditions  were 
spatially  variable,  the  central  place  patterns  would  then  be 
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distorted  to  account  for  that  variability.  Drawing  from  the 
work  of  Tobler,  Rushton  sought  to  develop  a  method  for 
deriving  expected  patterns  based  on  the  spatially  varying 
conditions.  Specif ically ,  he  wished  to  determine  the 
location  and  number  of  supply  points,  such  that  the  proximal 
area  equates  to  the  area  of  threshold  demand,  while 
simultaneously  minimizing  patron  travel  distance.  His 
approach  utilizes  a  distance  minimizing  concept  for  boundary 
definition.  He  suggests  several  possible  extensions  of  his 
method.  One  is  for  the  definition  of  tributary  regions  with 
population  and  distance  constraints.  For  example,  define 
the  regions  so  that  they  have  equal  populations  and  that 
no  individual  would  have  to  travel  further  than  some 
critical  distance  to  a  service  point.  He  calls  another 
extension  revised  boundary  definition.  In  cases  where  people 
have  free  choice  of  which  facility  to  patronize,  a 
space-preference  function  could  be  included  so  as  to  locate 
service  points  which  draw  equal  populations. 


CHAPTER  FIVE 


MODEL  SPECIFICATION 


No  existing  spatial  model  adequately  represents  the 


process  of  generating  composite  protection  zones.  Using 


the  practical  problem  of  protection  zone  delineation  and 


existing  spatial  models  generated  as  a  basis,  the 


assumptions  and  arguments  for  a  proposed  model  are  defined. 


A  fully  specified  model  is  incorporated  into  a  methodology 


for  solving  the  groundwater  protection  zone  problem.  It  is 


hoped  that  the  model  improves  the  method  for  definition  of 


composite  protection  zones,  as  well  as  yields  insights  into 


the  spatial  concepts  underlying  its  application. 


The  process  to  be  modeled  can  be  classified  as  an 


area-generating  point  process.  The  growth  model  concept  of 


Johnson  and  Mehl  (1939)  serves  as  the  foundation  for  the 


proposed  model.  This  model  generates  areal  patterns  or 


configurations  resulting  from  the  spread  of  influence  across 


an  area,  outward  from  a  set  of  points.  In  this  application, 


the  set  of  points  is  defined  by  well  locations  and  the 


composite  protection  zone  constitutes  the  areal  pattern. 


Assumptions  and  Arguments 


There  are  a  number  of  assumptions  which  will  be 


included  in  the  proposed  model  which  are  similar  to  those  of 


the  growth  model.  The  natural  representation  of  the 
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two-dimensional  growth  pattern  for  a  single  point  is  assumed 
to  be  a  circle  centered  at  that  same  point.  The  pattern  is 
produced  by  uniform  growth  in  all  directions  from  the  point. 
The  process  is  assumed  to  be  competitive  for  required 
resources  which  are  distributed  in  two-dimensional  space. 

The  environment  comprising  this  space  is  presumed  to  be 
uniform.  Unlike  the  growth  model,  an  assumption  is  included 
which  restricts  the  number  of  points  considered  to  two. 

This  supposition  is  expanded  to  accept  multiple  point 
processes  in  later  discussion. 

As  growth  about  a  single  point  results  exclusively  in 
circles,  the  arguments  specified  are  applicable  to  two-point 
patterns  whose  circles  contact  or  overlap.  Initially,  the 
circles  defining  the  space  required  by  each  point  are 
allowed  to  overlap.  Prom  this  configuration,  three  pieces 
of  information  concerning  the  process  are  extracted  for 
later  use  in  the  model.  The  total  area  of  the  final 
configuration  is  the  sum  of  the  areas  of  the  individual 
circles.  The  area  by  which  the  initial  figure  is 
"deficient"  is  equal  to  the  overlap  area  or  the  area  of 
intersection  between  the  two  circles.  The  initial 
configuration  boundary  is  used  as  the  figure  to  be  distorted 
in  the  growth  process. 

Boundary  Distortion 

Next,  a  condition  prohibiting  overlap  is  imposed  on  the 
process.  This  necessitates  expansion  of  the  initial 
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configuration  to  include  an  area  equivalent  in  size  to  the 
overlap  area.  The  distorted  figure  must  totally  enclose  the 
initial  figure.  Therefore,  if  the  configuration  boundary  is 
distorted,  it  must  be  extended  outward  from  the  interior  of 
the  figure. 

The  determination  of  the  boundary  distortion  technique 
is  the  focal  point  of  this  study.  Drawing  from  the 
previously  described  works  of  Rushton  (1972)  and  Tobler 
(1963),  a  unique  approach  is  posed.  Rushton  suggested  that 
a  space-preference  function  could  be  incorporated  into  map 
transformations  which  utilize  distance-minimizing  concepts 
for  boundary  definition.  The  space-preference  function  of 
interest  here  is  one  which  extends  the  boundary  nearer  the 
overlap  region  more  than  it  does  at  distances  farther  from 
the  overlap.  A  function  of  this  type  is  compatible  with  the 


inverse  distance  relationship  alluded  to  in  the  discussion 


of  the  study  objectives. 


The  distance  minimization  can  be  easily  understood  if 


one  views  the  expansion  process  as  occurring  about  a  single 


point  or  origin.  The  aggregate  distance  from  the  origin  to 


all  points  on  the  boundary  can  be  minimized.  Tobler  (1963) 


describes  map  transformations  which  are  analogous  to 


processes  of  this  type.  He  indicates  that  the  patterns 


generated  exhibit  radial  scale  distortion  from  the  origin. 


In  Tobler' s  work,  the  origin  is  the  point  from  which 
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migration  is  assumed  to  begin,  whereas  Rushton's  origin  was 
the  location  of  the  central  place  under  consideration. 

As  the  application  of  a  radial  scale  distortion 
transformation  to  the  study  model  requires  an  origin,  one 
must  be  selected.  However,  in  the  conceptualization  of  the 
present  model,  there  are  two  circle  center  points  surrounded 
by  the  overlap  configuration.  How  the  transformation  origin 
is  selected  has  a  definite  impact  on  the  transformed  figure. 

On  first  inspection,  one  might  assume  that  the  weighted 
mean  center  (using  circle  radius  as  weight)  of  the  two 
circle  centers  is  the  appropriate  distortion  origin.  But 
for  a  pair  of  circles  of  very  dissimilar  size,  the  weighted 
mean  center  lies  closer  to  the  center  of  the  larger  circle. 
Radial  distortion  from  this  point  would  not  reflect  the 
spatial  competition  pressure  present  in  the  overlap  area. 

The  distortion  would  be  almost  uniform  about  the  larger 
circle  center  with  a  slight  deviation  toward  the  small 
circle. 

Utilization  of  the  geographic  center  of  the  two  circle 
centers  as  the  transform  origin  provides  a  slight 
improvement  over  the  weighted  mean  center.  Nevertheless, 
dissimilar  circle  sizes  still  result  in  distortion  with 
exaggerated  influence  from  the  larger  circle.  For  this 
reason,  the  distortion  origin  selected  is  the  centroid 
(geometric  center)  of  the  overlap  region.  It  serves  to 
distort  the  boundary  ( relieve  competition  pressure)  relative 


to  the  distance  from  the  overlap.  Use  of  the  overlap 
centroid  then  seems  the  most  appropriate  alternative  for  the 
distortion  transformation. 


Using  the  overlap  centroid  for  the  transf ormation,  how 
is  the  distorted  figure  to  be  generated?  Given  that  the 
distorted  figure  must  totally  enclose  the  configuration,  a 
fundamental  concept  of  the  model  can  be  derived.  For  each 
point  on  the  transformed  boundary,  the  distance  to  the 
overlap  centroid  is  always  greater  than  or  equal  to  the 
colinear  distance  from  the  initial  boundary  to  the  same 
centroid  (Figure  7).  This  means  that  the  distortion 
transformation  may  add  to,  but  never  subtract  from  the 
boundary/centroid  interpoint  distance. 

Determining  the  distance  to  add  to  this  interpoint 
distance  is  the  key  to  modeling  this  spatial  process.  An 
inverse  distance  relation  has  been  discussed  several  times, 
but  never  specified.  The  interpoint  distance  (d)  could  be 
extended  by  a  factor  proportional  to  1/d  or  for  that  matter 
any  polynomial  function  in  d. 

Transformation  Function 

The  function  selected  for  use  in  this  study  is  one 
which  extends  interpoint  distance  by  a  factor  of  1/d2,  The 
function  was  not  chosen  arbitrarily,  but  selected  after 
contemplation  of  its  use  in  existing  models.  For  example, 
Newton's  Principle  of  Universal  Gravitation  and  Coulomb's 
Principle  of  Electrostatic  Attraction  both  calculate  forces 
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Figure  7:  Boundary/Centroid  Interpoint  Distance 
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as  inversely  proportional  to  the  squared  distance 
separating  the  objects  under  consideration.  Incorporating 
this  function  into  the  distortion  model  seemed  appropriate, 
especially  if  one  envisions  the  pressure  of  spatial 
competition  as  a  force.  When  this  pressure  is  diffused 
through  space  the  force  diminishes  in  such  a  manner  as  to  be 
inversely  proportional  to  the  square  of  the  distance  from 
the  original  location.  The  squared  term  arises  from  the 
increase  in  area  available  to  dissipate  the  force  as  one 
moves  away  from  a  point  of  higher  pressure.  This  increase 
is  analogous  to  the  effect  of  increasing  the  radius  of  a 
circle.  The  areal  increase  is  proportional  to  the  radial 
increase  squared.  The  point  from  which  the  diffusion 
process  takes  place  is  specified  as  the  overlap  centroid  for 
this  study. 

In  order  to  control  the  degree  of  distortion  introduced 
by  a  transformation,  a  scaling  factor  must  be  incorporated 
into  the  model.  This  factor  will  limit  the  maximum  distance 
extension  increment  to  a  desired  value.  Combining  the 
scaling  factor  and  the  distance  function  the  transformation 
function  takes  the  following  form: 

d'  =  d  +  K/d2 

where 

d'  =  transformed  distance  from  the  centroid  to  the 
boundary 


lit. 
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d  *  initial  distance  from  the  centroid  to  the 
boundary 

K  =  scaling  factor  to  limit  maximum  extension 
increment 

The  value  of  K  is  determined  by  inserting  the  shortest 
boundary/ interpoint  distance  for  a  given  configuration  into 
the  above  equation  as  d.  The  maximum  desired  increment  is 
added  to  d  to  get  a  value  for  d' .  K  can  then  be  determined 
algebraically. 

The  transformation  function  is  applied  to  each  point 
which  defines  the  initial  configuration  boundary.  This 
procedure  has  the  effect  of  expanding  the  configuration 
boundary.  The  areal  increase  attributable  to  the  distortion 
process  is  dependent  on  the  value  of  K,  as  well  as  the 
initial  overlap  configuration.  Since  a  single 
transformation  may  not  increase  the  figure  area  to  the  total 
area  required,  multiple  transformations  may  be  necessary. 

Using  an  iterative  process,  the  transformation  can  be 
repeated  until  the  distorted  figure  boundary  encompasses  an 
area  equal  to  (or  vary  within  a  selected  tolerance  of)  the 
total  area  requirement.  Each  repetition  requires 
recalculation  of  K  and  utilization  of  the  most  recently 
transformed  boundary  points.  The  distortion  process  is 
complete  when  the  area  requirement  is  met. 
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The  distortion  model  is  designed  to  accommodate  three 
possible  structural  overlap  configurations  (Figure  8). 
Lateral  overlap  is  the  simplest  case  and  involves  two 
circles  overlapping  side-to-side .  The  distortion  process 
described  to  this  point  considered  a  lateral  overlap 
structure  only. 

Embedded  overlap  occurs  when  one  circle  lies  entirely 
within  another.  For  this  two-circle  structure,  a  slight 
modification  of  the  distortion  model  is  necessary.  The 
centroid  of  overlap  is  now  defined  as  the  center  of  the 
embedded  circle.  The  initial  configuration  boundary  is  the 
outer  circle.  The  minimum  interpoint  distance  used  to 
determine  K  is  the  shortest  distance  from  the  centroid  to 
the  outer  circle. 

The  most  complex  structures  are  those  involving 
multiple  overlap.  They  consist  of  three  or  more  circles 
involving  any  combination  of  lateral  and/or  embedded 
overlap.  Model  modification  is  also  required  for  multiple 
overlap  cases.  Each  overlapping  pair  of  circles  is  treated 
as  a  separate  overlap  case  of  its  respective  simple  type 
(lateral  or  embedded).  The  outer  configuration  boundary 
(union)  of  the  multiple  overlap  is  used  for  all  scaling 
calculations  and  transformations.  The  distortion  process  is 
conducted  by  rotating  through  a  sequence  of  transformation 
origins,  one  for  each  overlap  pair.  A  single  transformation 
is  completed  for  one  origin,  then  the  areal  increase  is 
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compared  with  the  respective  origin's  overlap  area.  If  the 
area  requirement  is  met,  the  origin  is  dropped  from  the 
sequence.  Irrespective  of  whether  the  origin  is  dropped,  the 
transformation  origin  is  advanced  to  the  next  overlap 
centroid  in  sequence  and  the  transformation  process  is 
repeated.  The  procedure  is  continued  iteratively  until  all 
origins  are  dropped  from  the  sequence,  indicating  that  the 
total  area  requirement  has  been  met.  This  transform  - 
rotate  -  transform  -  rotate  sequence  allows  for  uniform 
influence  of  all  overlap  pairs  in  the  distortion  process. 


CHAPTER  SIX 


MODEL  APPLICATION 

The  goal  to  produce  an  objective  method  for  defining 
groundwater  protection  zones  is  of  primary  concern  in  this 
study.  The  manual  cartographic  procedure  previously 
described  falls  well  short  of  meeting  that  goal.  Kaplan  et 
al .  (1986)  indicate  that  a  great  deal  of  aquifer  mapping  is 

presently  done  by  hand  using  similar  labor-intensive, 
subjective  procedures.  They  suggest  that  utilizing  a 
geographic  information  system  (GIS)  would  "allow  faster, 
more  reliable,  and  scientifically  credible  results." 

Broten  et  al .  (1986)  also  support  the  use  of  a  GIS  for 
spatial  analyses  related  to  groundwater  contamination 
investigations.  They  cite  a  number  of  GIS  characteristics 
which  make  such  systems  desirable  for  groundwater  research. 
The  efficiency  and  accuracy  of  a  GIS  allows  for  simple 
automation  of  many  procedures.  Routine  functions  necessary 
for  spatial  modeling  are  also  included  in  these  systems. 

Data  management  and  map  preparation  are  greatly  enhanced  by 
the  use  of  a  GIS.  The  advantages  offered  by  these  systems 
are  not  strictly  confined  to  groundwater  research,  but 
rather  to  any  field  involving  spatial  analysis. 

To  capitalize  on  these  advantages,  a  GIS  was  used  in 
this  study.  The  software  employed  in  the  research  was  the 


ARC/INFO,  Revision  4.0  (Environmental  Systems  Research 
Institute)  geographic  information  system.  The  host  computer 
was  a  Prime  750  with  a  Fortran  77  compiler.  An  Altek  9000 
digitizing  table  was  used  as  the  data  input  device.  A 
Hewlett  Packard  7586B  pen  plotter  was  used  for  graphics 
output.  The  hardware  and  software  are  the  property  of  the 
Tallahassee  District  Office  of  the  U.S.  Geological  Survey, 
Tallahassee,  Florida. 

In  the  following  discussion  the  term  "system"  is  used 
to  refer  to  the  components  listed  above.  The  system  is  used 
in  any  way  possible  when  it  simplifies  or  enhances  the 
modeling  process.  While  describing  the  procedures  of  the 
process,  the  specific  system  commands  and  sequences  are 
omitted  and  only  the  general  descriptions  of  each  step  are 
discussed.  As  the  system  itself  does  not  have  all  the 
functional  capabilities  required  to  implement  the  model, 
additional  routines  had  to  be  developed.  These  routines  are 
Fortran  programs  written  to  perform  the  transformation 
process.  The  transformation  procedures  are  discussed  in 
greater  detail  than  the  system  functions  as  they  comprise 
the  essence  of  the  study. 

The  entire  modeling  process  consists  of  three  basic 
steps  to  produce  a  composite  groundwater  protection  zone. 
First,  utilizing  the  GIS,  input  data  for  the  transformation 
program  are  generated.  These  data  are  entered  into  the 
program  during  execution.  The  transformation  output  is  then 


used  as  input  into  the  system  for  generation  of  graphic 
output.  A  discussion  of  each  of  these  steps  follows.  First 
however,  a  description  of  the  transformation  algorithm  is 
presented  as  a  preview  to  understanding  the  required  inputs 

Transformation  Algorithm 

Rushton  (1972)  describes  the  use  of  locational 
coordinate  transformations  for  processes  requiring  revised 
boundary  definition.  He  suggests  that  by  imposing 
constraints  based  on  areal  size  and/or  space  preference 
functions,  a  transformation  can  be  used  to  define  a  new 
boundary  from  the  set  of  coordinates  describing  an  existing 
boundary.  The  solution  is  an  iterative  process  which  maps 
each  boundary  point  to  a  new  trial  point.  The  amount  and 
direction  of  each  move  is  dependent  on  the  constraints 
imposed.  This  type  of  transformation  procedure  was  used  as 
the  foundation  for  the  composite  protection  zone  generation 
process . 

In  Chapter  Five,  Model  Specification,  a  number  of  the 
conceptual  aspects  of  the  coordinate  transformation  are 
discussed.  The  transformation  origin  was  selected  to  be  the 
geometric  center  of  the  overlap  region.  The  space 
preference  function  chosen  was  one  which  transforms  each 
boundary  coordinate  to  a  new  point,  based  on  a  factor 
inversely  proportional  to  the  squared  distance  from  the 
original  boundary  point  to  the  overlap  centroid.  The  ereal 


constraint  imposed  was  the  sum  of  the  areas  of  the 
individual  protection  zones  involved. 

The  transformation  algorithm  was  derived  through 
algebraic  manipulation  of  the  general  linear  equation  and 
distance  formula.  Using  the  coordinates  of  the  centroid 
(Xc,Yc)  and  a  boundary  point  (X,Y)  to  be  transformed,  the 
equation  (slope  and  intercept)  of  a  line  through  the  two 
points  can  be  determined  in  the  form: 

Y  =  a  +  bX 

Using  the  distance  formula,  the  distance  (d)  between  the 
same  two  points  can  be  calculated  from: 

d  =  [(X-Xc)2  +  (y-Yc)2]* 

As  discussed  in  Chapter  Five,  the  transformation  serves  to 
extend  distance  d  to  a  new  distance  d'  by  the  function: 

d'  =  d  +  K/d2 

The  distance  each  point  is  extended  from  the  origin  is  then 
o 

equal  to  K/d  ,  with  K  being  the  previously  described  scaling 
factor . 

K  is  inserted  to  control  the  maximum  distance  extension 

possible  for  a  single  transformation.  Its  value  is 

2 

determined  by  setting  K/d  equal  to  the  maximum  extension 
desired  when  d  is  equal  to  the  shortest  initial 
centroid-boundary  distance  (see  Figure  7).  For  this  reason, 
the  shortest  distance  must  be  determined  prior  to  the 
transformation  and  be  included  as  input  to  the  process. 

This  value  remains  the  same  through  one  complete 
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transformation  iteration,  which  processes  each  coordinate 
pair  once.  The  maximum  extension  desired  for  the  protection 
zone  application  is  ten  feet.  This  figure  is  based  on  the 
U.S.  Geological  Survey  requirement  to  round  all  protection 
zone  radius  calculations  to  the  nearest  ten  feet.  The 
ten-foot  increment  also  seems  reasonable  in  order  that 
iterative  expansion  of  the  protection  zone  boundary  fall 
within  the  five-percent  areal  tolerance.  A  larger  increment 
might  cause  the  transformation  to  overshoot  the  areal 
requirement  beyond  the  tolerance  level. 

Once  the  extended  distance  d*  is  calculated  for  a 
boundary  point,  the  transformed  point  coordinates  (X1 ,Y' ) 
can  be  determined  by  the  following  derived  equations: 

X*  =  +  d’/tb2*!)*  +  xc 
Y’  =  a  +  b  X’ 

where 

d’  *  extended  distance  from  centroid  to 
boundary 

a  =*  intercept  of  line  through  centroid  and 
boundary  point 
b  =  slope  of  same  line 
Xc  =  X  coordinate  of  origin  (centroid) 

The  +  for  the  X1  term  arises  from  the  distance 
calculation  which  by  definition  results  in  only  positive 
values.  Since  some  boundary  coordinates  will  have  smaller  X 
values  for  the  boundary  points  than  for  the  origin  (negative 
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X  distance),  the  direction  of  distance  extension  must  be 
considered.  Thus,  for  all  boundary  points  with  an  X 
coordinate  value  less  than  Xc,  the  negative  value  of  d  is 
used  to  calculate  the  transformed  coordinate  X1 .  One  step 
in  the  derivation  of  the  transformed  coordinate  equations 
required  translation  of  the  centroid  coordinates  to  the 
zero,  zero  origin  (0,0).  The  addition  of  Xc  in  the  derived 
equation  for  X'  is  required  to  retranslate  the  values  back 
into  the  original  coordinate  system.  Neither  the  +  term, 
nor  the  retranslation  term  appear  in  the  equation  for  the 
transformed  Y'  value.  The  Y*  coordinate  is  determined  from 
the  general  linear  equation  using  X' .  Therefore,  the 
direction  of  the  distance  extension  and  origin  translation 
are  already  accounted  for. 

The  transformed  coordinates  X'  and  Y'  then  replace  the 
original  input  coordinates  X  and  Y.  The  transformation 
process  is  repeated  for  the  next  X,  Y  input  pair  until  all 
boundary  points  have  been  extended.  This  cycle  of 
processing  each  boundary  point  one  time  constitutes  one 
iteration  of  the  transformation.  After  one  iteration,  the 
value  of  the  shortest  distance  from  the  origin  to  the  new 
boundary  is  saved  for  further  iterations  if  necessary.  The 
new  value  is  simply  the  smallest  extended  distance  (d1) 
calculated  during  the  previous  iteration.  Additional 
iterations  are  required  if  the  areal  increase  due  to  the 
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transformation  is  less  than  the  overlap  area.  Then  overlap 
area  is  also  required  as  an  input. 

Preparation  of  Inputs 

Of  the  three  basic  steps  to  the  delineation  of  a 
composite  protection  zone,  the  input  preparation  is 
the  only  one  which  varies  and  is  dependent  on  the  overlap 
structure  under  consideration  (Figure  8).  Initially, 
preparation  of  the  data  for  a  two-well  lateral  overlap 
structure  is  described.  The  modifications  for  embedded  and 
multiple  structures  is  presented  later  in  this  section. 

The  source  for  all  inputs  to  the  transformation  program 
is  USGS  1:24,000  topographic  quadrangles  on  which  well 
locations  and  individual  protection  zones  are  plotted.  Also 
plotted  on  the  quadrangles  are  the  manually  created 
composite  protection  zones.  Ho  verification  of  either  well 
locations  or  protection  zone  radii  was  attempted,  since  a 
comparison  of  the  manual  and  transformed  composite  zones  was 
intended.  Correction  of  well  locations  and/or  individual 
zone  sizes  would  certainly  alter  the  resulting  transformed 
composite  protection  zone  and  therefore  make  comparisons 
meaningless . 

The  well  locations  are  digitized  into  the  system  and 
circular  polygons  are  generated  to  represent  the  individual 
protection  zones.  The  areas  of  the  two  circular  polygons  as 
determined  by  the  system  are  summed  to  provide  the  total 
areal  requirement  for  the  composite  protection  zone.  The 
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system  is  then  used  to  determine  the  intersection  of  the 
individual  protection  zones.  The  area  of  the  lens-shaped 
intersection  polygon  is  that  area  which  must  be  accounted 
for  by  the  transformation.  This  intersection  area  (denoted 
as  overlap)  is  another  input  required  by  the  transformation. 

The  G1S  is  also  used  to  determine  the  coordinates  of 
the  centroid  of  the  overlap  area  which  serves  as  the 
transformation  origin.  A  system  function  that  generates 
label  points  for  unlabeled  polygons  is  used  for  this 
purpose.  According  to  the  ARC/INFO  documentation,  the  label 
point  generated  for  each  polygon  would  be  located  at  its 
calculated  geometric  center.  Visual  inspection  of  the  label 
locations  for  the  lens-shaped  overlap  polygons  indicated 
that  the  algorithm  for  centroid  calculation  is  satisfactory 
for  these  particular  figures.  However,  additional  tests  of 
the  centroid  calculation  routine  for  irregular-shaped 
polygons  yielded  some  unpredictable  and  incorrect  geometric 
center  locations.  The  generated  centroids  are  used  for  the 
lens-shaped  polygons  since  they  are  the  only  shape  possible 
for  the  overlap  polygons  in  the  study.  The  system's 
centroid  determination  algorithm  may  be  inappropriate  if  the 
model  is  modified  to  include  overlap  polygons  which  are  not 
symmetrical . 

The  two  circular  polygons  representing  the  individual 
protection  zones  are  overlaid  using  the  system  to  obtain  the 
union  of  the  figures.  The  area  of  the  union  is  recorded  for 


input  into  the  transformation  as  the  initial  figure  area. 
Using  the  boundary  of  the  union  polygon  and  the  calculated 
overlap  centroid,  a  distance  measuring  function  is  invoked. 
Its  purpose  is  to  determine  the  shortest  distance  from  the 
centroid  to  the  boundary  of  the  polygon,  another  required 
input.  Finally,  a  system  function  is  utilized  to  generate  a 
file  consisting  of  X,Y  coordinates  which  define  the  union 
polygon's  boundary.  This  is  the  initial  set  of  coordinates 
transformed  in  the  composite  protection  zone  generation. 

At  this  point,  all  data  inputs  are  known,  but  they 
still  must  be  put  into  proper  input  format.  An  input 
parameter  file  is  created  which  includes:  the  X  and  Y 
coordinates  of  the  centroid,  the  shortest  distance  from  the 
centroid  to  the  union  boundary,  and  the  size  of  the  overlap 
area,  respectively.  One  final  parameter  included  in  this 
file  is  the  maximum  distance  extension  increment,  which  is 
previously  specified  as  ten  feet.  Since  all  coordinates  and 
distances  are  measured  in  digitizer  inches  on  a  1:24,000 
scale  map,  the  figure  for  the  maximum  increment  is  0.005 
Inches,  the  equivalent  of  ten  feet  on  the  earth's  surface. 

For  overlap  structures  other  than  the  lateral  type, 
modifications  of  the  input  preparation  procedure  are 
required.  Embedded  structures  actually  require  less  work 
for  gathering  inputs  than  does  the  lateral  overlap  case.  No 
intersection  process  is  necessary  as  the  embedded  circle  is 
the  intersection  polygon.  It  follows  that  the  area  and 


center  of  this  circle  are  the  area  and  centroid  of  the 
overlapping  polygon.  The  union  polygon  area  and  boundary 
are  simply  that  of  the  larger  individual  protection  zone  in 
the  embedded  structure.  All  other  inputs  are  obtained  in 
the  same  manner  described  for  the  lateral  overlap  case. 

The  multiple  overlap  structure  requires  additional 
input  data  and,  therefore,  additional  work  in 
pretransformation  data  processing.  The  individual 
protection  zones  are  considered  for  all  possible  pairs  of 
intersection.  For  purposes  of  the  study,  a  visual 
determination  of  overlap  occurrence  was  made.  If  a  large 
number  of  wells  is  involved,  a  6IS  analysis  may  be  required 
to  identify  all  existing  overlap  situations.  For  each 
overlap  condition  which  does  exist,  inputs  are  gathered  by 
the  methods  described  for  the  respective  simple  overlap 
structure  (lateral  or  embedded).  The  resulting  parameter 
file  consists  of  one  record  for  each  paired  overlap 
occurrence.  Each  record  contains  centroid  coordinates, 
shortest  distance  to  the  union  boundary,  overlap  area,  and 
the  maximum  extension  increment.  The  initial  figure  area, 
total  required  area,  and  union  boundary  are  determined  in  a 
manner  similar  to  that  previously  described  using  all  wells 
and  protection  zones  involved. 

Transformation  Program 
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software  consists  of  a  main  program  and  two  subroutines 
written  in  Fortran  77  code  (see  Appendix) .  An  overview  of 
each  of  the  three  program  segments  is  presented  in  the 
following  discussion. 

The  initial  step  of  the  main  program  requires 
interactive  input  of  the  parameter  file,  initial  figure 
area,  total  required  area,  boundary  coordinate  file,  and  an 
output  file  for  the  transformed  coordinates.  Each  variable 
for  each  record  in  the  parameter  file  is  read  into  an  array 
to  be  used  during  the  transformation  processing.  The  number 
of  centroids  input  (and  therefore  number  of  records)  is 
counted  and  displayed  as  a  check  on  the  input  process.  Each 
centroid  input  serves  as  transformation  origin  for 
coordinate  transformation.  The  overlap  area  associated  with 
each  centroid  is  summed  to  determine  the  amount  of  area 
which  must  be  added  to  the  union  figure. 

Adding  the  sum  of  the  overlap  area  and  the  initial 
figure  area,  the  result  is  equal  to  the  summed  areas  of  the 
individual  protection  zones  in  those  cases  where  only  two 
wells  overlap  the  same  area.  As  the  approach  used  for  the 
study  considered  pair-wise  overlap  only,  occurrence  of 
overlapping  areas  by  three  or  more  individual  protection 
zones  would  yield  an  artificially  high  value  for  the  sum  of 
overlap  areas  and  the  initial  figure  area.  Areas  of  three- 
zone  intersection  would  be  added  one-too-many  times,  four- 
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zone  intersection  two-too-many ,  and  so  forth,  in  overlap 
summation. 

If  the  transformation  was  executed  using  this  improper 
value,  the  resulting  composite  protection  zone  would  be 
larger  than  required.  For  this  reason  a  compensation  factor 
is  introduced  which  reduces  each  value  of  overlap  stored  in 
the  parameter  array  by  a  factor  of: 

(Total  Area  Required  -  Initial  Area) /Overlap  Sum 
This  correction  produces  a  transformed  figure  with  the 
desired  area  (sum  of  the  area  of  individual  protection 
zones )  . 

The  boundary  coordinates  for  the  union  figure  are  read 
into  another  array  and  the  number  of  vertices  is  counted  and 
displayed.  The  main  program  then  begins  a  loop  to  transform 
the  coordinates  using  the  first  record  (centroid)  in  the 
parameter  file.  A  subroutine  to  do  the  transformation  is 
called  and  executed  as  previously  described.  Next,  an 
area-calculating  subroutine  is  called  to  determine  and 
display  the  area  of  the  transformed  figure.  The  area  added 
to  the  previous  figure  is  calculated  and  subtracted  from  the 
overlap  area  associated  with  the  transformation  centroid  and 
stored  back  into  the  parameter  array.  This  allows  for 
flagging  a  centroid  when  its  overlap  area  has  been  accounted 
for  by  transformation  with  the  centroid  as  the  origin. 

The  procedure  continues  with  the  next  parameter  record 
in  sequence  (if  a  multiple  overlap  structure).  The  same 
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processes  are  repeated  as  the  transformed  boundary  is 
expanded  from  the  second  centroid.  After  a  single 
transformation  is  completed  for  each  centroid  in  the 
parameter  file,  one  iteration  of  the  process  is  complete. 
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The  processing  continues  from  the  first  parameter  record 
through  the  same  sequence.  Any  time  a  parameter  record  is 
read  and  its  overlap  value  is  zero,  the  centroid  is  skipped 
in  the  transformation  sequence.  When  all  the  overlap 
values  in  the  parameter  file  are  equal  to  zero,  the 
transformation  process  is  complete. 
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The  final  transformed  coordinates  contained  in  the 
coordinate  array  are  written  to  the  output  file.  The  number 
of  iterations  which  were  required  to  generate  the  composite 
protection  zone  is  displayed.  This  completes  the  execution 
of  the  transformation  program. 

Use  of  Output 

The  output  file  containing  the  transformed  coordinates 
requires  minor  format  modification  prior  to  importation  into 
the  GIS.  The  import  routine  is  followed  by  a  system 
function  invoked  to  build  polygon  structure.  This  function 


provides  an  area  calculation  for  the  composite  protection 
zone  which  is  used  as  a  check  for  the  transformation  program 
calculation.  Using  the  graphics  capabilities  of  the  system, 
the  composite  protection  zone  generated  is  viewed  on  the 
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CRT.  Hardcopy  output  is  then  produced  using  the  system 
plotter . 


► 


The  wells  selected  for  use  in  the  study  were  chosen 
from  those  requiring  composite  zones  presented  in  Vecchioli, 
et  al .  ,  (1987).  Those  selected  involved  a  combination  of 
the  three  possible  overlap  structures.  The  number  of  wells 
involved  in  any  single  composite  zone  was  restricted  to  five 
or  less  to  reduce  input  preparation  time.  Table  1  lists  the 
wells  used  in  the  study  along  with  associated  information 
about  each,  including:  the  quadrangle  on  which  the  well  is 
found,  well  ID,  latitude,  longitude,  individual  protection 
zone  radius,  and  the  generated  figure  number. 
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TABLE  1 


Wells  Selected  for  Composite 
Protection  Zone  Generation 
(Source:  Vecchioli  et  al.,  1987) 


Quad  Sheet 

Well 

Lat 

Lona 

Radius ( f t ) 

Fiqure 

Weekiwachee 

95 

28 

32 

27 

82 

34 

00 

1380 

9 

II 

1585 

28 

32 

39 

82 

34 

01 

1150 

9 

Crystal  River 

51 

28 

54 

06 

82 

34 

38 

1950 

11 

II 

52 

28 

54 

04 

82 

34 

38 

1510 

11 

Chassahowitzka 

#8 

28 

42 

23 

82 

31 

25 

1570 

13 

II 

#9 

28 

42 

12 

82 

31 

15 

1610 

13 

II 

#14 

28 

42 

30 

82 

31 

34 

1650 

13 

II 

#15 

28 

42 

34 

82 

31 

44 

1650 

13 

Citrus  Park 

1525 

28 

06 

12 

82 

31 

53 

1190 

15 

II 

1526 

28 

06 

19 

82 

31 

57 

1050 

15 

II 

1529 

28 

05 

45 

82 

32 

17 

990 

15 

II 

1530 

28 

05 

52 

82 

31 

58 

1360 

15 

II 

1531 

28 

05 

50 

82 

31 

57 

1330 

15 

Citrus  Park 

1537 

28 

03 

49 

82 

33 

17 

2240 

17 

II 

1538 

28 

03 

40 

82 

32 

33 

2250 

17 

II 

1539 

28 

03 

32 

82 

32 

40 

2250 

17 

II 

1540 

28 

03 

21 

82 

31 

54 

2260 

17 

Port  Richey 

434 

28 

16 

43 

82 

41 

49 

970 

19 

II 

435 

28 

16 

31 

82 

41 

51 

960 

19 

II 

436 

28 

16 

29 

82 

41 

53 

1140 

19 

II 

437 

28 

16 

27 

82 

41 

29 

1320 

19 
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CHAPTER  SEVEN 


DISCUSSION  OP  RESULTS 


The  results  of  the  study  are  discussed  in  the  two 


following  sections.  The  first  section  describes  the  program 


operation  characteristics.  The  computer-related  aspects  of 


the  composite  protection  zone  generation  procedure  with 


respect  to  data  processing  time  and  resource  requirements 


are  addressed.  A  discussion  of  the  generated  composite 


zones  with  an  analysis  of  each  figure  is  also  provided. 


Each  composite  protection  zone  generated  is  compared  with 


its  manually  delineated  equivalent  and  evaluated  with 


respect  to  the  defined  study  objectives. 


Program  Operation 


The  preparation  of  input  data  for  the  transformation 


program  accounted  for  a  substantial  majority  of  the  time 


required  to  generate  a  composite  protection  zone. 


Preparation  time  varied  from  one  to  four  hours  depending  on 


the  particular  situation.  The  time  increased  with  an 


increasing  number  of  wells,  but  depended  particularly  on  the 


number  of  paired  overlap  occurrences.  The  more  overlap 


involved,  the  greater  the  time  required  to  obtain  inputs. 


This  owes  to  the  fact  that  each  paired  overlap  occurrence 


required  determination  of  its  centroid  coordinates,  shortest 


distance  to  the  boundary,  and  overlap  area  for  program 


input . 
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The  size  of  the  boundary  file  used  as  input  for  the 
transformation  program  varied  according  to  the  well 
configuration  involved.  The  minimum  possible  number  of 
vertices  is  360.  A  simple  embedded  structure  with  a 
circular  union  figure  would  have  this  number  of  vertices. 
This  is  a  system-induced  minimum  based  on  the  number  of 
points  used  to  define  a  polygon  produced  by  using  a  circle¬ 
generating  function.  The  Crystal  River  configuration  is  an 
example  consisting  of  the  minimum  possible  vertex  count. 

The  largest  number  of  vertices  contained  in  a  study 
configuration  input  file  was  859  for  Citrus  Park  1.  It 
should  also  be  noted  that  the  vertex  count  of  the  composite 
prote  tion  zone  created  is  identical  to  the  number  input. 

The  transformation  merely  modifies  each  point  location:  no 
new  vertices  are  generated. 

Like  vertex  count,  the  number  of  transformation 
iterations  required  to  produce  the  composite  zones  varied 
according  to  well  configuration.  The  total  area  of  overlap 
was  not  in  itself  the  sole  factor  in  determining  the 
required  number  of  iterations.  Rather,  total  overlap 
relative  to  the  initial  figure  area  gave  a  better  indication 
of  how  many  transformation  cycles  were  necessary.  The 
transformation  for  the  Weekiwachee  Springs  composite  zone 
required  35  iterations,  the  fewest  of  any  of  the  study 
figures.  Needing  160  iterations,  the  Citrus  P;r.rk  1 


configuration  took  the  greatest  number  to  complete  the 
transformation. 

The  computer  time  required  to  generate  the  composite 
protection  zones  was  not  a  significant  factor  for  any  of  the 
study  configurations.  With  all  processing  times  of  less 
than  one  minute,  the  computer  resources  utilized  during 
transformation  were  negligible  in  comparison  with  those 
necessary  for  input  preparation.  This  vast  difference  in 
time  required  for  input  preparation  and  transformation 
processing  suggests  that  modifications  to  the  transformation 
program  might  be  feasible.  Program  changes  which 
incorporate  input  preparation  steps  into  the  transformation 
software  may  prove  beneficial  with  respect  to  time  savings. 
Such  modifications  are  discussed  in  the  following  chapter. 

Use  of  the  computer  for  modeling  the  spatial  allocation 
process  did  serve  to  meet  a  primary  objective  of  the 
research.  The  objective  of  developing  a  formal  method  for 
the  generation  of  unique  composite  protection  zones  was 
satisfied.  The  method's  ability  to  consistently  reproduce 
results  was  also  verified.  A  thorough  examination  of  output 
files  produced  by  separate  executions  of  the  program 
revealed  identical  coordinate  sets.  Additional  checks  were 
made  on  the  number  of  iterations  required  and  generated 
protection  zone  area  to  ensure  that  by  repeating  the 
procedure,  the  same  values  were  obtained. 
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Generated  Composite  Protection  Zones 

Each  of  the  configurations  used  for  composite 
protection  zone  generation  is  discussed  in  this  section. 

The  figures  depicting  the  generated  composite  zones  are 
included.  Also  depicted  on  the  figures  are  the  well 
locations  and  individual  overlapping  protection  zones  for 
which  the  composite  zone  was  produced.  All  the  figures  are 
1:24,000  scale  portrayals.  The  composite  zones  generated 
were  visually  compared  with  the  U.S.  Geological  Survey's 
manually  created  equivalents.  The  USGS  composite  protection 
zone  figures  presented  are  extracted  portions  of  figures 
found  in  Vecchioli  et  al .  (1987). 

Figure  9  represents  the  Weekiwachee  Springs  composite 
protection  zone.  This  case  was  selected  as  an  example  of  a 
simple  lateral  overlap  structure.  The  resulting  composite 
zone  compares  favorably  with  the  USGS  version  (Figure  10) 
and  exhibits  only  minor  differences.  The  USGS  delineation 
does  not  differentiate  the  degree  of  extension  about  the 
dissimilar-sized  individual  zones  as  the  transformed  figure 
does.  Based  on  the  specified  distance  function  in  the 
model,  the  transformed  protection  zone  appears  as  expected. 
The  function  could  be  modified  to  become  more  distance- 
sensitive  by  inserting  d3  or  d4  in  place  of  d2 .  Such  a 
change  would  produce  a  composite  protection  zone  similar  to 
Figure  9,  but  display  even  greater  differentiation  in 
boundary  extension  from  bottom  to  top  of  the  protection 
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zone.  For  the  Weekiwachee  Springs  example,  all  of  the  study 
objectives  as  defined  appear  to  have  been  met. 

The  Crystal  River  example  (Figure  11)  was  selected 
because  of  its  embedded  overlap  structure.  From  this 
particular  case,  an  interesting  though  not  readily  apparent 
aspect  of  the  spatial  model  is  demonstrated.  Due  to 
differential  extension  of  the  boundary  points  based  on  the 
distance  function,  the  transformed  distances  from  centroid 
to  ail  vertices  converge  to  the  same  value.  Ultimately,  if 
significant  overlap  exists,  successive  transformations 
produce  a  boundary  comprised  of  points  equally  distant  from 
the  centroid.  By  definition  then,  the  resulting  composite 
protection  zone  is  a  circle  centered  on  the  overlap 
centroid.  Convergence  to  a  circle  is  evident  in  the  Crystal 
River  example.  The  manual  (Figure  12)  and  computer- 
generated  protection  zones  are  almost  indistinguishable  for 
this  particular  situation.  The  resulting  figure  also  meets 
the  stated  research  objectives. 

Figure  13  depicts  a  multiple  overlap  structure  of  wells 
located  on  the  Chassahowitzka  map  sheet.  During  efforts  to 
generate  this  particular  composite  protection  zone,  two 
shortcomings  of  the  methodology  were  revealed.  The  problems 
encountered  concerned  the  shortest  distance  value  and  total 
area  requirement  used  for  multiple  overlap  transformations. 

A  shortest  distance  (SDIST)  value  is  determined  for 
each  centroid  in  the  multiple  overlap  structure  and  input 
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via  the  parameter  file.  On  the  initial  iteration  from  the 
first  centroid  in  the  parameter  array,  the  system- 
determined  shortest  distance  is  properly  used  for  the 
transformation.  An  updated  SDIST  is  then  calculated  with 
respect  to  the  new  boundary  and  substituted  into  the 
parameter  array.  When  the  program  advances  to  the  second 
centroid,  the  transformation  proceeds  using  the  system- 
determined  SDIST  value  obtained  prior  to  any  boundary 
extension.  This  results  in  a  smaller  than  desired  boundary 
extension  increment  for  all  coordinates,  since  the  scaling 
factor  (K)  is  calculated  with  an  erroneously  low  SDIST 
value . 

Even  though  the  shortest  distance  is  recalculated  after 
transformation  from  its  associated  centroid,  its  value  is 
not  updated  to  reflect  boundary  changes  effected  by  prior 
transformations.  This  is  true  for  all  coordinate  processing 
which  follows  the  initial  transformation  in  a  multiple 
overlap  situation.  The  net  result  is  a  less  than  desired 
boundary  extension  for  all  but  the  first  coordinate 
transformation . 

The  severity  of  this  methodological  error  is  difficult 
to  determine  with  respect  to  effect  on  the  resulting 
composite  protection  zone.  Analysis  of  the  problem  seemed 
to  suggest  that  the  flaw  was  not  overly  critical.  The 
worst-case  possibility  of  any  of  the  study  figures  is  one  in 
which  five  prior  transformations  extend  the  actual  shortest 


distance  by  a  total  of  50  feet  (10  foot  maximum  increment 
times  5  transformations) .  The  scaling  factor  calculated 
using  710  feet  (the  smallest  system-determined  SDIST) , 
rather  than  760  feet,  would  produce  a  maximum  extension 
increment  of  8.7  feet. 

This  worst-case  scenario  assumes  that  all  of  the 
transformations  had  the  same  centroid  location.  That  is  the 
only  way  possible  for  five  boundary  extensions  to  result  in 
a  50-foot  addition  to  the  shortest  distance  of  a  sixth 
point.  Although  a  single  centroid  location  is  a 
possibility,  it  seems  very  unlikely  that  such  a  situation 
might  exist.  Even  if  the  centroid  locations  are  slightly 
offset  in  a  more  realistic  situation,  the  effect  on  the 
composite  zone  is  believed  to  be  minor.  The  expansion  from 
several  closely  located  origins  results  in  a  figure  very 
similar  to  one  generated  based  on  transformation  from  any 
one  of  the  individual  centroids.  The  number  of  iterations 
will  be  increased  because  of  the  decreased  maximum  extension 
increment,  but  the  final  result  will  be  nearly  identical. 

As  the  centroid  locations  become  more  distant  from  one 
another,  the  influence  of  prior  transformations  has  less 
effect  on  the  shortest  distance  values.  Therefore,  based  on 
this  analysis,  the  problem  associated  with  the  shortest 
distance  values  was  not  deemed  to  be  a  significant 
shortcoming  of  the  study.  Corrective  actions  for  the 
problem  are  suggested  in  the  following  chapter. 
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The  second  problem,  which  surfaced  during  generation  of 
the  Chassahowitzka  composite  zone,  was  one  previously 
described  involving  overlapping  areas  common  to  more  than 
two  individual  protection  zones.  The  method  for  summation 
of  the  paired  overlap  area  with  the  union  area  resulted  in  a 
figure  larger  than  the  total  area  requirement.  It  was  this 
problem  that  necessitated  introduction  of  the  compensation 
factor  described  in  Chapter  Six. 

Once  both  obstacles  to  multiple  overlap  transformation 
were  either  dismissed  as  inconsequential  or  remedied,  the 
Chassahowitzka  figure  was  created.  The  composite  zone 
generated  was  slightly  more  bulbous  near  the  centroids  than 
U.S.  Geological  Survey  version  (Figure  14).  This  seems 
plausible  based  on  the  distance  function  incorporated  into 
the  model.  Even  with  that  variation,  the  two  composite 
zones  were  remarkably  similar.  The  generated  Chassahowitzka 
protection  zone  was  also  acceptable  with  respect  to  the 
study  objectives. 


Citrus  Park  1  was  another  multiple  overlap  structure 
selected  for  protection  zone  generation  (Figure  15).  In 
this  example,  the  well  configuration  is  more  complex  than 
that  in  the  previous  case.  Because  of  the  additional 
complexity,  the  manual  delineation  technique  would 
presumably  be  more  difficult.  Evidence  of  this  is 


demonstrated  by  the  significant  differences  found  in  the 
generated  and  manually  created  versions  of  the  composite 
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protection  zones.  The  Geological  Survey's  composite  zone 
(Figure  16)  seemed  to  be  constructed  to  yield  a  symmetrical 
figure.  The  transformed  protection  zone  exhibits  a  more 
even  distribution  of  overlap  about  the  entire  union  figure. 
The  tendency  toward  circular  convergence  is  also  evident  in 
Figure  15,  because  of  to  the  large  areas  of  paired  overlap 
involved.  Like  all  previously  discussed  figures,  the 
generated  Citrus  Park  1  composite  protection  zone  adheres  to 
the  standards  set  in  the  study  objectives. 

Figure  17,  Citrus  Park  2,  was  selected  as  another 
unique  multiple  overlap  example.  Although  the 
transformation-produced  protection  zone  again  displayed  more 
even  overlap  distribution  about  the  union  figure,  it  greatly 
resembled  its  USGS  counterpart  (Figure  18). 

Close  examination  of  Figure  17  reveals  the  most  serious 
flaw  identified  in  the  composite  protection  zone  generation 
process.  A  point  on  the  composite  zone  boundary  was 
transformed  to  a  location  inside  the  union  figure  of  the 
individual  zones.  This  occurrence  violates  the  fundamental 
study  objective  that  the  generated  protection  zone  totally 
enclose  all  individual  zones  from  which  it  is  created.  The 
flaw  does  not  appear  to  cause  drastic  problems  for  the 
Citrus  Park  2  configuration,  but  suggests  the  potential  for 
more  serious  consequences  in  other  situations. 

The  total  enclosure  problem  was  a  completely 
unanticipated  result  of  the  research.  This  unforeseen 


Figure  16:  USGS  Citj 


possibility  escaped  detection  as  a  result  of  the  model 
development  process.  The  model  initially  considered 
addressed  a  two-circle  overlap  condition.  Boundary 
extension  from  the  overlap  centroid  in  the  two-circle  case 
would  always  result  in  a  totally  enclosed  union  figure. 
This  fact  was  accepted  as  a  given  argument  in  the  basic 
model.  When  expanded  to  consider  multiple  overlap 
structures,  the  argument  was  applied  to  the  modified  model 
and  assumed  valid  for  all  configurations.  A  more  in-depth 
review  of  the  argument's  applicability  to  multiple 
structures  might  have  ascertained  its  lack  of  validity. 


The  problem  depicted  in  Figure  17  seems  to  have  been 


produced  by  transformation  from  the  centroid  located  at  the 
figure  center.  This  particular  centroid  has  a  large 
associated  overlap  area  requiring  many  transformations. 
Convergence  to  circular  boundary  expansion  is  demonstrated 
for  the  configuration.  The  numerous  transformations 
effectively  "force"  the  composite  zone  boundary  inside  of  an 
individual  protection  zone. 

The  Port  Richey  example  (Figure  19)  was  selected 
because  of  its  perceived  potential  for  producing  a  similar 
enclosure  problem.  The  configuration  did  contain  two 
relatively  large  overlap  areas  in  close  proximity  with  each 
other.  Numerous  expansions  from  the  two  centroids  could 
possibly  introduce  the  boundary  enclosure  problem. 
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Even  though  the  generated  composite  zone  was  very 
similar  to  the  U.S.  Geological  Survey's  version  (Figure  20), 
the  transformation  did  indeed  locate  composite  boundary 
points  inside  an  individual  protection  zone.  The  distance 
within  the  individual  zone  was  small,  but  it  resulted  in  a 
figure  which  fails  to  meet  all  of  the  study  objectives.  In 
addition,  the  problems  associated  with  the  last  two 
generated  figures  discussed  raise  concerns  about  the 
applicability  of  the  methodology  to  a  number  of  multiple 
overlap  structures. 


II 


CHAPTER  EIGHT 


CONCLUSIONS 

The  conclusions  drawn  from  this  study  are  presented  in 
the  two  following  sections.  The  first  section  is  a 
discussion  of  the  model's  applicability  to  the  problem  of 
groundwater  protection  as  well  as  other  potential  uses. 
Possible  modifications  of  the  methodology  are  also 
introduced.  The  second  section  provides  a  summary  of  the 
conclusions . 

Applicability  of  the  Model 

For  composite  protection  zone  generation,  the  modeling 
process  produced  results  which  satisfied  all  of  the 
objectives  in  four  of  the  six  examples.  The  simple  cases 
involving  only  two  wells  with  either  lateral  or  embedded 
structures,  provided  pleasing  results  when  compared  with 
their  manually  created  equivalents.  In  those  cases  where 
the  figures  differed,  the  computer-generated  form  is  deemed 
to  be  a  more  reasonable  representation  of  the  composite 
zone,  if  only  because  of  the  formality  of  the  generation 
procedure.  It  is  believed  that  for  any  two-well  overlap 
structure,  the  transformed  composite  zone  is  preferable  to 
its  manually  delineated  equivalent. 

The  boundary  enclosure  problem  described  earlier  clouds 
the  issue  of  the  model's  applicability  for  multiple  overlap 
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structures.  In  two  of  the  four  multiple  overlap  examples, 
the  resulting  composite  zone  was  totally  acceptable  in  terms 
of  the  study  objectives  and  figure  appearance.  Two  other 
multiple  overlap  structures  produced  composite  zones  which 
did  not  totally  enclose  all  individual  protection  zones 
involved  in  the  configuration.  Thus,  the  model's 
applicability  requires  verification  on  a  case-by-case  basis. 
Since  an  unspecified  intention  of  the  study  was  to  develop  a 
universally  applicable  model,  additional  research  involving 
the  model's  modification  is  suggested. 

Other  potential  areas  exist  for  application  of  model. 

In  any  situation  where  objects  are  in  competition  for  space, 
the  model  may  have  utility.  Whether  the  spatial  competition 
process  involves  biological,  environmental,  or  economic 
phenomena,  the  proposed  model  may  provide  an  acceptable 
research  tool.  As  with  the  groundwater  applications, 
consideration  should  be  given  to  model  modification  when  the 
process  concerns  more  than  two  objects. 

Possible  Modifications 


A  number  of  modifications  to  the  composite  protection 
zone  generation  process  are  possible.  The  changes  discussed 
are  ones  which  would  alleviate  those  problems  which  are 
considered  to  be  merely  annoyances  as  well  as  major  flaws. 
Several  alterations  in  the  transformation  program  would 
result  in  extensive  reduction  in  the  time  required  to 
prepare  the  data  inputs.  Insertion  of  a  subroutine  to 


calculate  the  intersection  (overlap)  area  of  the  individual 
protection  zones  would  eliminate  a  tedious  portion  of  input 
preparation.  The  algorithm  for  the  subroutine  is  not  simple 
(Rosewell,  1986),  but  could  easily  be  incorporated  into  the 
program. 

A  centroid  determination  subroutine  might  also  be 
included  in  the  program.  Any  one  of  a  number  of  existing 
algorithms  could  be  used.  The  result  would  be  removal  of 
the  GIS  centroid  determination  step  from  the  input 
preparation  procedures. 

Another  modification  to  the  program  concerns  the 
shortest  distance  (SDIST)  determination.  As  mentioned  in 
Chapter  Seven,  the  value  of  SDIST  used  during  multiple 
overlap  transformation  is  correct  only  for  the  first 
transformation.  Since  the  subroutine  for  shortest  distance 
calculation  already  exists  in  the  transformation  software, 
the  problem  can  easily  be  corrected.  Changing  the  order  of 
the  program's  SDIST  calculation  to  a  point  immediately  prior 
to  transformation  processing,  rather  than  following  it, 
would  eliminate  the  problem.  The  reordering  of  processing 
steps  also  reduces  input  preparation  time,  as  a 
system-determined  shortest  distance  value  would  no  longer  be 
required  as  input. 

Each  of  the  changes  suggested  to  this  point  would  have 
desirable  effects  on  input  preparation  time  and  program 
operation.  In  fact,  it  is  estimated  that  the  maximum  input 
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data  preparation  time  required  for  any  of  the  configurations 
addressed  in  this  study,  could  be  reduced  from  four  hours  to 
less  than  one.  At  the  same  time,  program  execution  is  not 
expected  to  exceed  three  minutes,  a  two-minute  increase  over 
the  unmodified  execution  time.  The  modifications  do 
decrease  reliance  on  GIS-supplied  capabilities ,  since  they 
would  be  written  into  the  program.  This  contradicts  the 
previously  stated  intention  to  utilize  routine  spatial 
functions  available  within  the  system.  However,  the 
cumbersome,  time-consuming  nature  of  GIS  function  use 
dictates  that  the  modifications  be  incorporated  into  the 
program  and  that  system  reliance  be  reduced. 

Another  modification  should  be  included  to  enhance  the 
program  integrity  concerning  the  compensation  factor. 

Recall  that  the  factor  is  applied  to  each  paired  overlap 
area  when  a  configuration  contains  three  or  more  individual 
protection  zones  which  intersect  the  same  area.  It  serves 
as  a  correction  to  the  total  area  requirement  determined  by 
adding  the  summed  overlap  areas  tc  the  area  of  the  union. 

It  is  the  intersection  areas  common  to  more  than  two 
individual  zones  which  produce  the  summation  problem. 
Therefore,  areas  comprised  exclusively  of  two-zone  overlap 
should  not  be  reduced  by  the  compensation  factor. 

It  is  unclear  as  to  which  is  the  best  method  of 
handling  the  problem.  One  possibility  is  to  consider  each 
overlap  occurrence  based  on  the  number  of  individual  zones 
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which  intersect  to  form  an  overlap  area.  This  would 
eliminate  the  lens-shaped  polygon  of  paired  overlap  as  the 
only  possibility  considered.  The  multiple  intersection 
polygons  could  have  very  irregular  shapes.  As  noted, 
experience  with  the  system  indicates  that  the  GIS  centroid 
calculation  routine  may  not  prove  satisfactory  for  such 
configurations.  Since  an  improved  centroid  calculation 
routine  was  suggested  for  an  earlier  modification,  the 
system's  shortcoming  is  not  viewed  as  a  problem. 

The  result  of  this  change  is  that  each  polygon  formed 
by  overlaying  all  individual  zones  is  treated  as  a  separate 
transformation  entity.  Knowledge  of  the  location  of  its 
centroid,  its  area,  and  the  number  of  individual  zones  which 
intersect  to  form  it  is  required  for  transformation.  The 
area  added  from  a  particular  centroid  during  transformation 
would  be  based  on  the  number  of  zones  which  intersected  to 
form  the  associated  overlap  polygon.  An  area  equivalent  to 
that  of  the  polygon  is  added  for  two-zone  intersection.  The 
intersection  of  three,  four,  or  five  individual  protection 
zones  requires  addition  of  two,  three,  and  four  times  the 
intersection  polygon  area  respectively.  Comparison  of 
figures  generated  using  a  program  incorporating  this 
modification  with  those  produced  by  the  existing 
transformation  would  provide  for  interesting  research. 

One  additional  modeling  error  also  mandates  program 
modification.  The  boundary  enclosure  problem,  considered  to 


be  a  major  program  flaw,  necessitates  refinement  of  the 
model  to  accommodate  all  possible  structures. 

Unfortunately,  an  acceptable  approach  to  solving  the 
enclosure  dilemma  has  not  been  devised.  Extended  research 
into  the  problem  is  required.  It  is  believed  that  with  such 
efforts,  the  program  could  potentially  be  made  applicable  to 
any  and  all  overlap  structures. 

Summary 

The  spatial  model  and  its  associated  transformation 
program  appear  to  provide  an  acceptable  approach  to  defining 
composite  groundwater  protection  zones  for  simpler  overlap 
structures.  In  addition,  application  of  the  model  to 
analogous  spatial  processes  is  entirely  feasible.  Two- 
object  competition  for  space  can  be  modeled  using  the 
program  if  the  competitive  process  adheres  to  the 
assumptions  and  arguments  specified. 

A  number  of  modifications  are  also  possible  for  the 
program.  Several  of  the  changes  would  greatly  reduce  the 
time  required  to  prepare  the  inputs  for  the  transformation 
process.  Other  program  alterations  are  recommended  to 
increase  transformation  integrity  and  applicability.  The 
changes  in  the  shortest  distance  determination  and 
compensation  factor  application  would  result  in  a  program 
which  is  more  consistent  with  study  objectives. 

The  most  severe  flaw  associated  with  the  model  is  the 


boundary  enclosure  problem.  The  problem  affected  the 


applicability  of  the  model  to  some  multiple  overlap 
structures.  Although  program  modifications  are  suggested, 
the  specific  changes  required  to  produce  universal 
applicability  in  multiple  overlap  situations  are  absent. 
Additional  research  concerning  the  model  and  transformation 
process  is  necessary  before  specific  changes  can  be  cited. 
However,  acceptability  of  the  existing  program  for  use  on 
multiple  structures  can  be  evaluated  on  an  individual  basis 
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APPENDIX  A 


***  PROGRAM  TO  GENERATE  COMPOSITE  PROTECTION  ZONES  *** 

♦★★VARIABLE  DECLARATION  *** 

CHARACTER  PARMF*20 , INBND*20 , 0UTBND*20 
INTEGER  NVERT, NEXT ( 50 ) ,J, ITER 

DOUBLE  PRECISION  XCENT(50) ,YCENT(50) ,SDIST(50) ,0VRLAP(50) 
DOUBLE  PRECISION  TAREA , AREA , SUBTRT , XIN( 2000 ) , YIN ( 2000 ) 
DOUBLE  PRECISION  ADDIST(50) , RQAREA , OVRSUM 

***  INTERACTIVE  INPUTS  *** 

WRITE (1,*)  'ENTER  PARAMETER  FILE' 

READ(1 , 30)  PARMF 

WRITE(1,*)  'ENTER  INITIAL  FIGURE  AREA' 

READ(1 , 25)  TAREA 

WRITE (1,*)  'ENTER  TOTAL  AREA  REQUIRED' 

READ(1 , 25)  RQAREA 
WRITE(1 ,*)  'ENTER  BOUNDARY  FILE' 

READ(1 , 30)  INBND 
WRITE (1,*)  'ENTER  OUTPUT  FILE' 

READ(1 , 30)  OUTBND 

***  OPENS  INPUT  FILES  *** 

OPEN ( 8 , FILE- PARMF , STATUS- ' OLD ' .ERR-999) 

OPEN ( 9 . FI LE-INBND , STATUS- ' OLD ' , ERR-9 9  9 ) 

OPEN (10, FILE -OUTBND , STATUS- ' NEW ' , ERR- 7 7  7 ) 

***  READ  PARAMETER  FILE  INTO  ARRAY  *** 

DO  10  K  -  1,50 

READ  (8, 70, END-15)  XCENT(K) ,YCENT(K) , SDIST(K) ,OVRLAP(K) , 
$  ADDIST(K) 

NEXT(K)  -  1 

I  CONTINUE 

i  NCENT  -  K  -  1 

WRITE  (1,*)  'NCENT  -’.NCENT 

C***  LOOP  TO  SUM  OVERLAP  AREAS  ♦** 

OVRSUM  -  0.0 
DO  90  K  -  1, NCENT 

OVRSUM  -  OVRSUM  +  OVRLAP(K) 

90  CONTINUE 

C 


$ 


c  *** 


75 

c 

c  *** 

c 


40 

c 

c  *** 

45 

c 

M  - 1— « _ «- 

L,  XXX 
*** 
c 

100 
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MULTIPLE  OVERLAP  COMPENSATION  FACTOR  *** 

DO  75  K  -  l.NCENT 

OVRLAP(K)  -  OVRLAP(K)  *(RQAREA  -  TAREA)  /  OVRSUM 
CONTINUE 

READ  INPUT  COORDS  INTO  ARRAY  *** 

DO  40  I  -  1,2000 

READ  (9,20, END-45)  XIN(I) ,YIN(I) 

CONTINUE 

CALCULATE  NUMBER  OF  COORD  PAIRS  *** 

NVERT  -  I  -  1 

WRITE  (1,*)  'EXIT  COORDS  ARRAY,  NVERT  -'.NVERT 

LOOP  TO  TRANSFORM  FROM  EACH  CENTROID  WITH  OVERLAP  REMAINING 


C 

C  *** 

c 

c 

c  *** 

c 

c 

c  *** 


c 

c  *** 

c 


ITER  -  0 

DO  55  J  -  l.NCENT 

IF  (NEXT(J)  .EQ.  0)  THEN 
GOTO  55 
END  IF 

ITER  -  ITER  +  1 

CALL  COORDINATE  TRANSFORMATION  SUBROUTINE  *** 

CALL  TRANS (ADDIST, SDIST, NVERT ,XIN, YIN, XCENT , YCENT , J ) 

CALL  AREA  CALCULATING  SUBROUTINE  *** 

CALL  CALCAR (XIN, YIN, NVERT, AREA) 

UPDATE  OVERLAP  REMAINING  AFTER  TRANSFORMATION  *** 

SUBTRT  -  AREA  -  TAREA 
WRITE  (1,*)  'AREA  -'.AREA 
OVRLAP(J)  -  OVRLAP(J)  -  SUBTRT 
TAREA  -  AREA 

SET  FLAG  TO  SKIP  CENTROID  WHEN  OVERLAP  ACCOUNTED  FOR  *** 

IF  (OVRLAP(J)  .LE.  0)  THEN 
NEXT(J)  -  0 

WRITE  (1,*)  'OVERLAP' ,J, ZERO' 

ENDIF 
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55  CONTINUE 
C 

c  ***  LOOP  TO  CHECK  IF  ALL  OVERLAP  IS  ACCOUNTED  FOR  AT  ALL 

CENTROIDS 

C 

DO  60  L  -  I.NCENT 

IF  (NEXT(L)  .EQ.  1)  THEN 
GOTO  100 
END  IF 

60  CONTINUE 

C 

C  ***  WRITE  FINAL  TRANSFORMED  COORDS  FROM  ARRAY  TO  OUTPUT  FILE  *** 
DO  80  I  -  l.NVERT 

WRITE  (10,20)  XIN(I) ,YIN(I) 

80  CONTINUE 

C 

WRITE  (1,*)  'NUMBER  OF  ITERATIONS ITER 
C 

C  ***  FORMAT  FOR  PARAMETER  FILE  *** 

70  FORMAT  (5(2X,F6.3)) 

C 

C  ***  FORMAT  FOR  INITIAL  FIGURE  AREA  *** 

25  FORMAT  (F6.3) 

C 

C  ***  FORMAT  FOR  INPUT  AND  OUTPUT  COORDS  FILES  *** 

20  FORMAT  ( 2 ( 6X , F9 . 6 ) ) 

C 

C  ***  FORMAT  FOR  INTERACTIVE  INPUT  FILE  NAMES  *** 

30  FORMAT  (A20) 

C 

C  ***  CLOSE  ALL  FILES  *** 

CLOSE  (8) 

CLOSE  (9) 

CLOSE  (10) 

C 

GOTO  1000 
C 

C  ***  FILE  OPENING  MESSAGES  *** 

999  WRITE  (6,*)  'ERROR  IN  OPENING  INPUT  FILES' 

777  WRITE  (6,*)  'ERROR  IN  OPENING  OUTPUT  FILE' 

1000  STOP 
END 

C 

C 

C  ***  SUBROUTINE  TO  TRANSFORM  COORDINATES  *** 

C 


KFACT 


C 

C  *** 

c 

c 

c  *** 
c 

c 

c  *** 
c 

c 

C  ilnHf 

c 

c  *** 
c 

c  *** 
c 

Q  *** 

c 

c  *** 


c 

c  *** 


c 

c  *** 
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SUBROUTINE  TRANS  ( ADDIST , SDIST , NVERT , XIN , YIN , XCENT , YCENT , J ) 
DOUBLE  PRECISION 

,  SHORT , SLOPE , INTERC , A, B , NEWDI S , XOUT , YOUT 
DOUBLE  PRECISION  ADDIST(50) , SDIST(50) ,XCENT(50) ,YCENT(50) 
DOUBLE  PRECISION  XIN(2000) ,YIN(2000) 

INTEGER  NVERT, J 

CALCULATES  EXPANSION  FACTOR  *** 

KFACT-ADDIST(J)*SDIST(J)**2 

INITIALIZE  SHORT  TO  A  LARGE  VALUE  *** 

SHORT— 100 . 0 

EXECUTE  COORDINATE  TRANSFORMATION  LOOP  *** 

DO  50  M-l, NVERT 

CALCULATE  SLOPE  OF  LINE  BETWEEN  INPUT  COORD  AND  CENTROID  *** 
SLOPE— (YIN(M) -YCENT(J) )/(XIN(M) -XCENT(J) ) 

CALCULATE  INTERCEPT  OF  LINE  *** 

INTERC-YCENT ( J ) - SLOPE*XCENT ( J ) 

CALCULATE  DISTANCE**2  BETWEEN  INPUT  COORD  AND  CENTROID  *** 
A— (XIN(M) -XCENT(J) )**2+(YIN(M) -YCENT(J) )**2 

CALCULATE  EXTENDED  DISTANCE  *** 

NEWDI S-SQRT (A ) +KFACT/A 

SAVES  SMALLEST  NEWDIST  INTO  SHORT  TO  UPDATE  PARAMETERS  *** 

IF  ( NEWDI S  .LT.  SHORT)  THEN 
SHORT— NEWD I S 
ENDIF 

MAKES  NEWDIS  NEGATIVE  FOR  XIN  VALUE  TO  LEFT  OF  XCENT  *** 

IF  (XCENT(J)  .GT.  XIN(M))  THEN 
NEWDIS— 1 . 0*NEWDIS 
ENDIF 

CALCULATES  DERIVED  TRANSFORM  FACTOR  *** 

B  -  SL0PE**2  +1.0 


$ 

'.‘1 

$ 


K; 


rvj v-av'.vv 
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C  ***  CALCULATES  OUTPUT  COORDS  FROM  TRANSFORMATION  *** 
XOUT-NEVDIS  /  SQRT(B)+XCENT(J) 

YOUT-INTERC+SLOPE*XOUT 

C 

C  ***  WRITE  OUTPUT  COORDS  TO  COORD  ARRAY  *** 

XIN(M)  -  XOUT 
YIN(M)  -  YOUT 
C 

50  CONTINUE 

C 

C  ***  SETS  SDIST  TO  SMALLEST  NEWDIST  AND  UPDATES  PARAMETER  ARRAY 
*** 

25  SDIST(J)— SHORT 

C 

RETURN 

END 

C 

C  ***  SUBROUTINE  TO  CALCULATE  AREA  *** 

C 

SUBROUTINE  CALCAR  (XIN, YIN, NVERT, AREA) 

C 

C  ***  VARIABLE  DECLARATION  AND  INITIALIZATION  *** 

DOUBLE  PRECISION  XIN (2000) ,YIN(2000) .AREA 
REAL  XOLD , YOLD , YORIG , X , Y 
AREA— 0 . 0 
XOLD— XIN(NVERT) 

YORIG-YIN( NVERT) 

YOLD— 0 . 0 
C 

C  ***  LOOP  TO  REPEAT  FOR  EACH  COORD  PAIR  *** 

DO  60  N-l, NVERT 
X-XIN(N) 

Y-YIN(N) -YORIG 

AREA-AREA+ ( XOLD - X) * ( YOLD+Y) 

XOLD— X 
YOLD-Y 
60  CONTINUE 
C 

AREA-ABS (0 . 5* AREA) 

RETURN 

END 


* 


